Author: [log in to unmask]
Date: Mon Jun 22 19:15:18 2015
New Revision: 3185
Log:
more cleanup
Modified:
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/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/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 22 19:15:18 2015
@@ -2,17 +2,13 @@
import hep.aida.IHistogram1D;
import hep.aida.IHistogram2D;
-import java.util.HashMap;
import java.util.List;
import java.util.Map;
import org.apache.commons.math.stat.StatUtils;
-import org.hps.recon.ecal.triggerbank.AbstractIntData;
-import org.hps.recon.ecal.triggerbank.TIData;
-import org.hps.recon.ecal.triggerbank.TestRunTriggerData;
+import org.hps.recon.ecal.cluster.ClusterUtilities;
import org.lcsim.event.CalorimeterHit;
import org.lcsim.event.Cluster;
import org.lcsim.event.EventHeader;
-import org.lcsim.event.GenericObject;
import org.lcsim.geometry.Detector;
/**
@@ -41,6 +37,7 @@
IHistogram1D clusterSizePlot;
IHistogram1D clusterEnergyPlot;
IHistogram1D clusterTimes;
+ IHistogram1D clusterTimeMean;
IHistogram1D clusterTimeSigma;
//mg...12/13/2014
IHistogram1D twoclusterTotEnergy;
@@ -55,8 +52,11 @@
IHistogram2D fiducialenergyVsX;
IHistogram2D energyVsY;
IHistogram2D energyVsX;
+ IHistogram2D energyVsT;
IHistogram2D xVsY;
IHistogram2D pairsE1vsE2;
+ IHistogram2D pairsT1vsT2;
+ IHistogram1D pairsDeltaT;
int nEvents = 0;
int nTotHits = 0;
@@ -66,7 +66,6 @@
double sumClusterEnergy = 0;
double sumClusterTime = 0;
boolean fillHitPlots = true;
- private Map<String, Double> monitoredQuantityMap = new HashMap<>();
String[] ecalQuantNames = {"avg_N_hits", "avg_Hit_Energy",
"avg_N_clusters", "avg_N_hitsPerCluster", "avg_Cluster_Energy", "avg_ClusterTime"};
double maxE = 1.1;
@@ -90,6 +89,7 @@
this.fillHitPlots = fill;
}
+ @Override
protected void detectorChanged(Detector detector) {
System.out.println("EcalMonitoring::detectorChanged Setting up the plotter");
aida.tree().cd("/");
@@ -105,14 +105,18 @@
clusterCountPlot = aida.histogram1D(plotClustersDir + triggerType + "/"+clusterCollectionName + " Cluster Count per Event", 10, -0.5, 9.5);
clusterSizePlot = aida.histogram1D(plotClustersDir + triggerType + "/"+clusterCollectionName + " Cluster Size", 10, -0.5, 9.5);
clusterEnergyPlot = aida.histogram1D(plotClustersDir + triggerType + "/"+clusterCollectionName + " Cluster Energy", 100, -0.1, maxE);
- clusterTimes = aida.histogram1D(plotClustersDir + triggerType + "/"+clusterCollectionName + " Cluster Time Mean", 200, 0, 4.0 * 50);
+ clusterTimes = aida.histogram1D(plotClustersDir + triggerType + "/"+clusterCollectionName + " Cluster Seed Times", 400, 0, 4.0 * 50);
+ clusterTimeMean = aida.histogram1D(plotClustersDir + triggerType + "/"+clusterCollectionName + " Cluster Time Mean", 400, 0, 4.0 * 50);
clusterTimeSigma = aida.histogram1D(plotClustersDir + triggerType + "/"+clusterCollectionName + " Cluster Time Sigma", 100, 0, 10);
twoclusterTotEnergy = aida.histogram1D(plotClustersDir + triggerType + "/"+clusterCollectionName + " Two Cluster Energy Sum", 100, 0, maxE);
twoclusterEnergyAsymmetry = aida.histogram1D(plotClustersDir + triggerType + "/"+clusterCollectionName + " Two Cluster Energy Asymmetry", 100, 0, 1.0);
- xVsY = aida.histogram2D(plotClustersDir + triggerType + "/"+clusterCollectionName + "X vs Y (NHits >1)", 200, -200.0, 200.0, 85, -85.0, 85.0);
+ energyVsT = aida.histogram2D(plotClustersDir + triggerType + "/"+clusterCollectionName + " Energy vs time", 400, 0.0, 200.0, 100, -0.1, maxE);
+ xVsY = aida.histogram2D(plotClustersDir + triggerType + "/"+clusterCollectionName + " X vs Y (NHits >1)", 200, -200.0, 200.0, 85, -85.0, 85.0);
energyVsX = aida.histogram2D(plotClustersDir + triggerType + "/"+clusterCollectionName + " Energy vs X", 50, 0, 1.6, 50, .0, 200.0);
energyVsY = aida.histogram2D(plotClustersDir + triggerType + "/"+clusterCollectionName + " Energy vs Y", 50, 0, 1.6, 50, 20.0, 85.0);
pairsE1vsE2 = aida.histogram2D(plotClustersDir + triggerType + "/"+clusterCollectionName + "Pair E1 vs E2", 50, 0, 2, 50, 0, 2);
+ pairsT1vsT2 = aida.histogram2D(plotClustersDir + triggerType + "/"+clusterCollectionName + "Pair T1 vs T2", 200, 0, 100, 200, 0, 100);
+ pairsDeltaT = aida.histogram1D(plotClustersDir + triggerType + "/" + clusterCollectionName + " Pair Time Difference", 100, -20.0, 20.0);
fiducialClusterCountPlot = aida.histogram1D(plotClustersDir + triggerType + "/"+plotFidCutDir + clusterCollectionName + " Cluster Count with Fiducal Cut", 10, -0.5, 9.5);
fiducialClusterSizePlot = aida.histogram1D(plotClustersDir+ triggerType + "/" +plotFidCutDir + clusterCollectionName + " Cluster Size with Fiducal Cut", 10, -0.5, 9.5);
@@ -159,10 +163,6 @@
List<Cluster> clusters;
if (event.hasCollection(Cluster.class, clusterCollectionName))
clusters = event.get(Cluster.class, clusterCollectionName);
- else if (event.hasCollection(Cluster.class, clusterCollectionName))
- clusters = event.get(Cluster.class, clusterCollectionName);
- else if (event.hasCollection(Cluster.class, clusterCollectionName))
- clusters = event.get(Cluster.class, clusterCollectionName);
else {
clusterCountPlot.fill(0);
return;
@@ -173,9 +173,11 @@
int fidcnt = 0;
for (Cluster cluster : clusters) {
clusterEnergyPlot.fill(cluster.getEnergy());
+ clusterTimes.fill(ClusterUtilities.findSeedHit(cluster).getTime());
+ energyVsT.fill(ClusterUtilities.findSeedHit(cluster).getTime(), cluster.getEnergy());
sumClusterEnergy += cluster.getEnergy();
double[] times = new double[cluster.getCalorimeterHits().size()];
- double[] energies = new double[cluster.getCalorimeterHits().size()];
+// double[] energies = new double[cluster.getCalorimeterHits().size()];
CalorimeterHit seed = cluster.getCalorimeterHits().get(0);
int ix = seed.getIdentifierFieldValue("ix");
int iy = seed.getIdentifierFieldValue("iy");
@@ -196,11 +198,11 @@
int size = 0;
for (CalorimeterHit hit : cluster.getCalorimeterHits()) {
- energies[size] = hit.getCorrectedEnergy();
+// energies[size] = hit.getCorrectedEnergy();
times[size] = hit.getTime();
size++;
}
- clusterTimes.fill(StatUtils.mean(times, 0, size));
+ clusterTimeMean.fill(StatUtils.mean(times, 0, size));
clusterSizePlot.fill(size); //The number of "hits" in a "cluster"
clusterTimeSigma.fill(Math.sqrt(StatUtils.variance(times, 0, size)));
sumHitPerCluster += size;
@@ -212,16 +214,18 @@
if (clusters.size() == 2) {
Cluster cl1 = clusters.get(0);
Cluster cl2 = clusters.get(1);
- double[] p1 = cl1.getPosition();
- double[] p2 = cl2.getPosition();
+// double[] p1 = cl1.getPosition();
+// double[] p2 = cl2.getPosition();
double e1 = cl1.getEnergy();
double e2 = cl2.getEnergy();
+ double t1 = ClusterUtilities.findSeedHit(cl1).getTime();
+ double t2 = ClusterUtilities.findSeedHit(cl2).getTime();
twoclusterTotEnergy.fill(e1 + e2);
twoclusterEnergyAsymmetry.fill(Math.abs(e1 - e2) / (e1 + e2));
pairsE1vsE2.fill(e1, e2);
-
- }
-
+ pairsT1vsT2.fill(t1, t2);
+ pairsDeltaT.fill(t1 - t2);
+ }
}
@Override
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 22 19:15:18 2015
@@ -84,10 +84,11 @@
IHistogram1D eneOverp;
IHistogram1D deltaXAtCal;
IHistogram1D deltaYAtCal;
- //IHistogram2D trackXvsECalX = aida.histogram2D(plotDir + finalStateParticlesColName + "/" + triggerType+ "/" + "track X vs ECal X", 50, -300, 300.0, 50, -300, 300.0);
- //IHistogram2D trackYvsECalY = aida.histogram2D(plotDir + finalStateParticlesColName + "/" + triggerType+ "/" + "track Y vs ECal Y", 50, -100, 100.0, 50, -100, 100.0);
+// IHistogram2D trackXvsECalX;
+// IHistogram2D trackYvsECalY;
IHistogram2D trackPvsECalE;
IHistogram2D trackTvsECalT;
+ IHistogram1D timeMatchDeltaT;
/* number of unassocaited tracks/event */
IHistogram1D nUnAssTracksHisto;
@@ -127,10 +128,11 @@
eneOverp = aida.histogram1D(plotDir + finalStateParticlesColName + "/" + triggerType + "/" + "Cluster Energy Over TrackMomentum", 50, 0, 2.0);
deltaXAtCal = aida.histogram1D(plotDir + finalStateParticlesColName + "/" + triggerType + "/" + "delta X @ ECal (mm)", 50, -50, 50.0);
deltaYAtCal = aida.histogram1D(plotDir + finalStateParticlesColName + "/" + triggerType + "/" + "delta Y @ ECal (mm)", 50, -50, 50.0);
- //IHistogram2D trackXvsECalX = aida.histogram2D(plotDir + finalStateParticlesColName + "/" + triggerType+ "/" + "track X vs ECal X", 50, -300, 300.0, 50, -300, 300.0);
- //IHistogram2D trackYvsECalY = aida.histogram2D(plotDir + finalStateParticlesColName + "/" + triggerType+ "/" + "track Y vs ECal Y", 50, -100, 100.0, 50, -100, 100.0);
+// trackXvsECalX = aida.histogram2D(plotDir + finalStateParticlesColName + "/" + triggerType + "/" + "track X vs ECal X", 50, -300, 300.0, 50, -300, 300.0);
+// trackYvsECalY = aida.histogram2D(plotDir + finalStateParticlesColName + "/" + triggerType + "/" + "track Y vs ECal Y", 50, -100, 100.0, 50, -100, 100.0);
trackPvsECalE = aida.histogram2D(plotDir + finalStateParticlesColName + "/" + triggerType + "/" + "track mom vs ECal E", 50, 0.1, beamEnergy * maxFactor, 50, 0.1, beamEnergy * maxFactor);
trackTvsECalT = aida.histogram2D(plotDir + finalStateParticlesColName + "/" + triggerType + "/" + "track T vs ECal T", 200, 0.0, 200.0, 100, -25.0, 25.0);
+ timeMatchDeltaT = aida.histogram1D(plotDir + finalStateParticlesColName + "/" + triggerType + "/" + "ECal T minus track T", 200, -25, 175);
/* number of unassocaited tracks/event */
nUnAssTracksHisto = aida.histogram1D(plotDir + finalStateParticlesColName + "/" + triggerType + "/" + "Number of unassociated tracks per event", 5, 0, 5);
}
@@ -256,10 +258,11 @@
deltaXAtCal.fill(dx);
deltaYAtCal.fill(dy);
/* here are some plots for debugging track-cluster matching */
- // aida.histogram2D(plotDir + finalStateParticlesColName + "/" + triggerType+ "/" + "track X vs ECal X").fill(trackPosAtEcal.x(), clusterPosition.x());
- // aida.histogram2D(plotDir + finalStateParticlesColName + "/" + triggerType+ "/" + "track Y vs ECal Y").fill(trackPosAtEcal.y(), clusterPosition.y());
+// trackXvsECalX.fill(trackPosAtEcal.x(), clusterPosition.x());
+// trackYvsECalY.fill(trackPosAtEcal.y(), clusterPosition.y());
trackPvsECalE.fill(fsPart.getMomentum().magnitude(), fsPart.getEnergy());
trackTvsECalT.fill(ClusterUtilities.getSeedHitTime(fsCluster), TrackUtils.getTrackTime(fsTrack, hitToStrips, hitToRotated));
+ timeMatchDeltaT.fill(ClusterUtilities.getSeedHitTime(fsCluster) - TrackUtils.getTrackTime(fsTrack, hitToStrips, hitToRotated));
// if(dy<-20)
// System.out.println("Big deltaY...")
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 22 19:15:18 2015
@@ -9,17 +9,14 @@
import hep.physics.vec.Hep3Vector;
import hep.physics.vec.VecOp;
import java.util.ArrayList;
-import java.util.Collection;
import java.util.List;
import java.util.Map.Entry;
+import org.hps.recon.tracking.TrackUtils;
import org.lcsim.event.EventHeader;
-import org.lcsim.event.LCRelation;
import org.lcsim.event.ReconstructedParticle;
import org.lcsim.event.RelationalTable;
import org.lcsim.event.Track;
-import org.lcsim.event.TrackerHit;
import org.lcsim.event.Vertex;
-import org.lcsim.event.base.BaseRelationalTable;
import org.lcsim.geometry.Detector;
/**
@@ -31,8 +28,6 @@
*/
public class TridentMonitoring extends DataQualityMonitor {
- private final String helicalTrackHitRelationsCollectionName = "HelicalTrackHitRelations";
- private final String rotatedHelicalTrackHitRelationsCollectionName = "RotatedHelicalTrackHitRelations";
private double ebeam = 1.05;
String finalStateParticlesColName = "FinalStateParticles";
String unconstrainedV0CandidatesColName = "UnconstrainedV0Candidates";
@@ -41,8 +36,7 @@
String trackListName = "MatchedTracks";
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"};
- boolean debug = false;
- private String plotDir = "TridentMonitoring/";
+ private final String plotDir = "TridentMonitoring/";
IHistogram2D trackTime2D;
IHistogram1D trackTimeDiff;
IHistogram2D vertexMassMomentum;
@@ -95,6 +89,10 @@
int nPassTrkCuts = 0;
int nPassClusterCuts = 0;
+
+ public void setEbeam(double ebeam) {
+ this.ebeam = ebeam;
+ }
@Override
protected void detectorChanged(Detector detector) {
@@ -136,46 +134,47 @@
@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(Track.class, trackListName))
- return;
+ 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(Track.class, trackListName)) {
+ return;
+ }
//check to see if this event is from the correct trigger (or "all");
- if (!matchTrigger(event))
- return;
+ if (!matchTrigger(event)) {
+ return;
+ }
nRecoEvents++;
- RelationalTable hittostrip = new BaseRelationalTable(RelationalTable.Mode.MANY_TO_MANY, RelationalTable.Weighting.UNWEIGHTED);
- List<LCRelation> hitrelations = event.get(LCRelation.class, helicalTrackHitRelationsCollectionName);
- for (LCRelation relation : hitrelations)
- if (relation != null && relation.getFrom() != null && relation.getTo() != null)
- hittostrip.add(relation.getFrom(), relation.getTo());
-
- RelationalTable hittorotated = new BaseRelationalTable(RelationalTable.Mode.ONE_TO_ONE, RelationalTable.Weighting.UNWEIGHTED);
- List<LCRelation> rotaterelations = event.get(LCRelation.class, rotatedHelicalTrackHitRelationsCollectionName);
- for (LCRelation relation : rotaterelations)
- if (relation != null && relation.getFrom() != null && relation.getTo() != null)
- hittorotated.add(relation.getFrom(), relation.getTo());
+ RelationalTable hitToStrips = TrackUtils.getHitToStripsTable(event);
+ RelationalTable hitToRotated = TrackUtils.getHitToRotatedTable(event);
List<Track> trks = event.get(Track.class, trackListName);
int ntracks = trks.size();
- if (ntracks > nTrkMax || ntracks < 2)
- return;
+ if (ntracks > nTrkMax || ntracks < 2) {
+ return;
+ }
List<ReconstructedParticle> fspList = event.get(ReconstructedParticle.class, finalStateParticlesColName);
int npos = 0;
- for (ReconstructedParticle fsp : fspList)
- if (fsp.getCharge() > 0)
+ for (ReconstructedParticle fsp : fspList) {
+ if (fsp.getCharge() > 0) {
npos++;
- if (npos < 1 || npos > nPosMax)
- return;
+ }
+ }
+ if (npos < 1 || npos > nPosMax) {
+ return;
+ }
nPassBasicCuts++;//passed some basic event-level cuts...
@@ -185,47 +184,49 @@
// v0 & vertex-quality cuts
// Hep3Vector v0Mom = uncV0.getMomentum();
Hep3Vector v0Mom = VecOp.add(uncV0.getParticles().get(1).getMomentum(), uncV0.getParticles().get(0).getMomentum());
- if (v0Mom.z() > v0PzMax || v0Mom.z() < v0PzMin)
- break;
- if (Math.abs(v0Mom.y()) > v0PyMax)
- break;
- if (Math.abs(v0Mom.x()) > v0PxMax)
- break;
+ if (v0Mom.z() > v0PzMax || v0Mom.z() < v0PzMin) {
+ break;
+ }
+ if (Math.abs(v0Mom.y()) > v0PyMax) {
+ break;
+ }
+ if (Math.abs(v0Mom.x()) > v0PxMax) {
+ break;
+ }
Hep3Vector v0Vtx = uncVert.getPosition();
- if (Math.abs(v0Vtx.z()) > v0VzMax)
- break;
- if (Math.abs(v0Vtx.y()) > v0VyMax)
- break;
- if (Math.abs(v0Vtx.x()) > v0VxMax)
- break;
+ if (Math.abs(v0Vtx.z()) > v0VzMax) {
+ break;
+ }
+ if (Math.abs(v0Vtx.y()) > v0VyMax) {
+ break;
+ }
+ if (Math.abs(v0Vtx.x()) > v0VxMax) {
+ break;
+ }
List<Track> tracks = new ArrayList<Track>();
ReconstructedParticle electron = null, positron = null;
- for (ReconstructedParticle particle : uncV0.getParticles())
-// tracks.addAll(particle.getTracks()); //add add electron first, then positron...down below
- if (particle.getCharge() > 0)
+ for (ReconstructedParticle particle : uncV0.getParticles()) // tracks.addAll(particle.getTracks()); //add add electron first, then positron...down below
+ {
+ if (particle.getCharge() > 0) {
positron = particle;
- else if (particle.getCharge() < 0)
+ } else if (particle.getCharge() < 0) {
electron = particle;
- else
+ } else {
throw new RuntimeException("expected only electron and positron in vertex, got something with charge 0");
+ }
+ }
+ if (electron == null || positron == null) {
+ throw new RuntimeException("vertex needs e+ and e- but is missing one or both");
+ }
tracks.add(electron.getTracks().get(0));
tracks.add(positron.getTracks().get(0));
- if (tracks.size() != 2)
+ if (tracks.size() != 2) {
throw new RuntimeException("expected two tracks in vertex, got " + tracks.size());
+ }
List<Double> trackTimes = new ArrayList<Double>();
for (Track track : tracks) {
- int nStrips = 0;
- double meanTime = 0;
- for (TrackerHit hit : track.getTrackerHits()) {
- Collection<TrackerHit> htsList = hittostrip.allFrom(hittorotated.from(hit));
- for (TrackerHit hts : htsList) {
- nStrips++;
- meanTime += hts.getTime();
- }
- }
- meanTime /= nStrips;
- trackTimes.add(meanTime);
+ trackTimes.add(TrackUtils.getTrackTime(track, hitToStrips, hitToRotated));
}
trackTime2D.fill(trackTimes.get(0), trackTimes.get(1));
trackTimeDiff.fill(trackTimes.get(0) - trackTimes.get(1));
@@ -258,8 +259,9 @@
@Override
public void printDQMData() {
System.out.println("TridentMonitoring::printDQMData");
- for (Entry<String, Double> entry : monitoredQuantityMap.entrySet())
+ for (Entry<String, Double> entry : monitoredQuantityMap.entrySet()) {
System.out.println(entry.getKey() + " = " + entry.getValue());
+ }
System.out.println("*******************************");
}
@@ -278,8 +280,9 @@
@Override
public void printDQMStrings() {
for (int i = 0; i < 9; i++)//TODO: do this in a smarter way...loop over the map
-
+ {
System.out.println("ALTER TABLE dqm ADD " + fpQuantNames[i] + " double;");
+ }
}
IFitResult fitVertexPosition(IHistogram1D h1d, IFitter fitter, double[] init, String range) {
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 22 19:15:18 2015
@@ -17,21 +17,21 @@
import java.util.Map.Entry;
import java.util.logging.Level;
import java.util.logging.Logger;
-import org.hps.recon.particle.HpsReconParticleDriver;
import org.hps.recon.tracking.TrackUtils;
import org.hps.recon.vertexing.BilliorTrack;
import org.hps.recon.vertexing.BilliorVertex;
import org.hps.recon.vertexing.BilliorVertexer;
import org.lcsim.event.EventHeader;
import org.lcsim.event.ReconstructedParticle;
+import org.lcsim.event.RelationalTable;
import org.lcsim.event.Track;
import org.lcsim.event.Vertex;
import org.lcsim.geometry.Detector;
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
*
@@ -54,6 +54,29 @@
double sumVz = 0.0;
double sumChi2 = 0.0;
+ /* V0 Quantities */
+ /* Mass, vertex, chi^2 of fit */
+ /* unconstrained */
+ IHistogram1D unconMass;
+ IHistogram1D unconVx;
+ IHistogram1D unconVy;
+ IHistogram1D unconVz;
+ IHistogram1D unconChi2;
+ /* beamspot constrained */
+
+ IHistogram1D nV0;
+ IHistogram1D bsconMass;
+ IHistogram1D bsconVx;
+ IHistogram1D bsconVy;
+ IHistogram1D bsconVz;
+ IHistogram1D bsconChi2;
+ /* target constrained */
+ IHistogram1D tarconMass;
+ IHistogram1D tarconVx;
+ IHistogram1D tarconVy;
+ IHistogram1D tarconVz;
+ IHistogram1D tarconChi2;
+
IHistogram2D pEleVspPos;
IHistogram2D pEleVspPosWithCut;
IHistogram2D pyEleVspyPos;
@@ -90,8 +113,7 @@
IHistogram1D sumChargeHisto;
IHistogram1D numChargeHisto;
- boolean debug = false;
- private String plotDir = "V0Monitoring/";
+ private final String plotDir = "V0Monitoring/";
double beamEnergy = 1.05; //GeV
double maxFactor = 1.25;
@@ -114,25 +136,25 @@
/* V0 Quantities */
/* Mass, vertex, chi^2 of fit */
/* unconstrained */
- IHistogram1D unconMass = aida.histogram1D(plotDir + triggerType + "/" + "Unconstrained Mass (GeV)", 100, 0, 0.200);
- IHistogram1D unconVx = aida.histogram1D(plotDir + triggerType + "/" + "Unconstrained Vx (mm)", 50, -10, 10);
- IHistogram1D unconVy = aida.histogram1D(plotDir + triggerType + "/" + "Unconstrained Vy (mm)", 50, -10, 10);
- IHistogram1D unconVz = aida.histogram1D(plotDir + triggerType + "/" + "Unconstrained Vz (mm)", 50, -50, 50);
- IHistogram1D unconChi2 = aida.histogram1D(plotDir + triggerType + "/" + "Unconstrained Chi2", 25, 0, 25);
+ unconMass = aida.histogram1D(plotDir + triggerType + "/" + "Unconstrained Mass (GeV)", 100, 0, 0.200);
+ unconVx = aida.histogram1D(plotDir + triggerType + "/" + "Unconstrained Vx (mm)", 50, -10, 10);
+ unconVy = aida.histogram1D(plotDir + triggerType + "/" + "Unconstrained Vy (mm)", 50, -10, 10);
+ unconVz = aida.histogram1D(plotDir + triggerType + "/" + "Unconstrained Vz (mm)", 50, -50, 50);
+ unconChi2 = aida.histogram1D(plotDir + triggerType + "/" + "Unconstrained Chi2", 25, 0, 25);
/* beamspot constrained */
- IHistogram1D nV0 = aida.histogram1D(plotDir + triggerType + "/" + "Number of V0 per event", 10, 0, 10);
- IHistogram1D bsconMass = aida.histogram1D(plotDir + triggerType + "/" + "BS Constrained Mass (GeV)", 100, 0, 0.200);
- IHistogram1D bsconVx = aida.histogram1D(plotDir + triggerType + "/" + "BS Constrained Vx (mm)", 50, -10, 10);
- IHistogram1D bsconVy = aida.histogram1D(plotDir + triggerType + "/" + "BS Constrained Vy (mm)", 50, -10, 10);
- IHistogram1D bsconVz = aida.histogram1D(plotDir + triggerType + "/" + "BS Constrained Vz (mm)", 50, -50, 50);
- IHistogram1D bsconChi2 = aida.histogram1D(plotDir + triggerType + "/" + "BS Constrained Chi2", 25, 0, 25);
+ nV0 = aida.histogram1D(plotDir + triggerType + "/" + "Number of V0 per event", 10, 0, 10);
+ bsconMass = aida.histogram1D(plotDir + triggerType + "/" + "BS Constrained Mass (GeV)", 100, 0, 0.200);
+ bsconVx = aida.histogram1D(plotDir + triggerType + "/" + "BS Constrained Vx (mm)", 50, -10, 10);
+ bsconVy = aida.histogram1D(plotDir + triggerType + "/" + "BS Constrained Vy (mm)", 50, -10, 10);
+ bsconVz = aida.histogram1D(plotDir + triggerType + "/" + "BS Constrained Vz (mm)", 50, -50, 50);
+ bsconChi2 = aida.histogram1D(plotDir + triggerType + "/" + "BS Constrained Chi2", 25, 0, 25);
/* target constrained */
- IHistogram1D tarconMass = aida.histogram1D(plotDir + triggerType + "/" + "Target Constrained Mass (GeV)", 100, 0, 0.200);
- IHistogram1D tarconVx = aida.histogram1D(plotDir + triggerType + "/" + "Target Constrained Vx (mm)", 50, -1, 1);
- IHistogram1D tarconVy = aida.histogram1D(plotDir + triggerType + "/" + "Target Constrained Vy (mm)", 50, -1, 1);
- IHistogram1D tarconVz = aida.histogram1D(plotDir + triggerType + "/" + "Target Constrained Vz (mm)", 50, -10, 10);
- IHistogram1D tarconChi2 = aida.histogram1D(plotDir + triggerType + "/" + "Target Constrained Chi2", 25, 0, 25);
+ tarconMass = aida.histogram1D(plotDir + triggerType + "/" + "Target Constrained Mass (GeV)", 100, 0, 0.200);
+ tarconVx = aida.histogram1D(plotDir + triggerType + "/" + "Target Constrained Vx (mm)", 50, -1, 1);
+ tarconVy = aida.histogram1D(plotDir + triggerType + "/" + "Target Constrained Vy (mm)", 50, -1, 1);
+ tarconVz = aida.histogram1D(plotDir + triggerType + "/" + "Target Constrained Vz (mm)", 50, -10, 10);
+ tarconChi2 = aida.histogram1D(plotDir + triggerType + "/" + "Target Constrained Chi2", 25, 0, 25);
pEleVspPos = aida.histogram2D(plotDir + triggerType + "/" + "P(e) vs P(p)", 50, 0, beamEnergy * maxFactor, 50, 0, beamEnergy * maxFactor);
pEleVspPosWithCut = aida.histogram2D(plotDir + triggerType + "/" + "P(e) vs P(p): Radiative", 50, 0, beamEnergy * maxFactor, 50, 0, beamEnergy * maxFactor);
pyEleVspyPos = aida.histogram2D(plotDir + triggerType + "/" + "Py(e) vs Py(p)", 50, -0.04, 0.04, 50, -0.04, 0.04);
@@ -171,31 +193,39 @@
@Override
public void process(EventHeader event) {
/* make sure everything is there */
- if (!event.hasCollection(ReconstructedParticle.class, finalStateParticlesColName))
+ if (!event.hasCollection(ReconstructedParticle.class, finalStateParticlesColName)) {
return;
- if (!event.hasCollection(ReconstructedParticle.class, unconstrainedV0CandidatesColName))
+ }
+ if (!event.hasCollection(ReconstructedParticle.class, unconstrainedV0CandidatesColName)) {
return;
- if (!event.hasCollection(ReconstructedParticle.class, beamConV0CandidatesColName))
+ }
+ if (!event.hasCollection(ReconstructedParticle.class, beamConV0CandidatesColName)) {
return;
- if (!event.hasCollection(ReconstructedParticle.class, targetV0ConCandidatesColName))
+ }
+ if (!event.hasCollection(ReconstructedParticle.class, targetV0ConCandidatesColName)) {
return;
+ }
//check to see if this event is from the correct trigger (or "all");
- if (!matchTrigger(event))
+ 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) {
Vertex uncVert = uncV0.getStartVertex();
- aida.histogram1D(plotDir + triggerType + "/" + "Unconstrained Vx (mm)").fill(uncVert.getPosition().x());
- aida.histogram1D(plotDir + triggerType + "/" + "Unconstrained Vy (mm)").fill(uncVert.getPosition().y());
- aida.histogram1D(plotDir + triggerType + "/" + "Unconstrained Vz (mm)").fill(uncVert.getPosition().z());
- aida.histogram1D(plotDir + triggerType + "/" + "Unconstrained Mass (GeV)").fill(uncV0.getMass());
- aida.histogram1D(plotDir + triggerType + "/" + "Unconstrained Chi2").fill(uncVert.getChi2());
-
- aida.histogram2D(plotDir + triggerType + "/" + "Mass vs Vz").fill(uncV0.getMass(), uncVert.getPosition().z());
+ unconVx.fill(uncVert.getPosition().x());
+ unconVy.fill(uncVert.getPosition().y());
+ unconVz.fill(uncVert.getPosition().z());
+ unconMass.fill(uncV0.getMass());
+ unconChi2.fill(uncVert.getChi2());
+
+ massVsVtxZ.fill(uncV0.getMass(), uncVert.getPosition().z());
VtxXVsVtxZ.fill(uncVert.getPosition().x(), uncVert.getPosition().z());
VtxYVsVtxZ.fill(uncVert.getPosition().y(), uncVert.getPosition().z());
VtxXVsVtxY.fill(uncVert.getPosition().x(), uncVert.getPosition().y());
@@ -222,23 +252,25 @@
}
double pe = ele.getMomentum().magnitude();
double pp = pos.getMomentum().magnitude();
- aida.histogram2D(plotDir + triggerType + "/" + "P(e) vs P(p)").fill(pe, pp);
- aida.histogram2D(plotDir + triggerType + "/" + "Px(e) vs Px(p)").fill(ele.getMomentum().x(), pos.getMomentum().x());
- aida.histogram2D(plotDir + triggerType + "/" + "Py(e) vs Py(p)").fill(ele.getMomentum().y(), pos.getMomentum().y());
+ pEleVspPos.fill(pe, pp);
+ pxEleVspxPos.fill(ele.getMomentum().x(), pos.getMomentum().x());
+ pyEleVspyPos.fill(ele.getMomentum().y(), pos.getMomentum().y());
if (pe < v0MaxPCut && pp < v0MaxPCut && (pe + pp) > v0ESumMinCut && (pe + pp) < v0ESumMaxCut)//enrich radiative-like events
- aida.histogram2D(plotDir + triggerType + "/" + "P(e) vs P(p): Radiative").fill(pe, pp);
+ {
+ pEleVspPosWithCut.fill(pe, pp);
+ }
}
List<ReconstructedParticle> beamConstrainedV0List = event.get(ReconstructedParticle.class, beamConV0CandidatesColName);
- aida.histogram1D(plotDir + triggerType + "/" + "Number of V0 per event").fill(beamConstrainedV0List.size());
+ nV0.fill(beamConstrainedV0List.size());
for (ReconstructedParticle bsV0 : beamConstrainedV0List) {
nTotV0++;
Vertex bsVert = bsV0.getStartVertex();
- aida.histogram1D(plotDir + triggerType + "/" + "BS Constrained Vx (mm)").fill(bsVert.getPosition().x());
- aida.histogram1D(plotDir + triggerType + "/" + "BS Constrained Vy (mm)").fill(bsVert.getPosition().y());
- aida.histogram1D(plotDir + triggerType + "/" + "BS Constrained Vz (mm)").fill(bsVert.getPosition().z());
- aida.histogram1D(plotDir + triggerType + "/" + "BS Constrained Mass (GeV)").fill(bsV0.getMass());
- aida.histogram1D(plotDir + triggerType + "/" + "BS Constrained Chi2").fill(bsVert.getChi2());
+ bsconVx.fill(bsVert.getPosition().x());
+ bsconVy.fill(bsVert.getPosition().y());
+ bsconVz.fill(bsVert.getPosition().z());
+ bsconMass.fill(bsV0.getMass());
+ bsconChi2.fill(bsVert.getChi2());
sumMass += bsV0.getMass();
sumVx += bsVert.getPosition().x();
sumVy += bsVert.getPosition().y();
@@ -249,32 +281,36 @@
List<ReconstructedParticle> targetConstrainedV0List = event.get(ReconstructedParticle.class, targetV0ConCandidatesColName);
for (ReconstructedParticle tarV0 : targetConstrainedV0List) {
Vertex tarVert = tarV0.getStartVertex();
- aida.histogram1D(plotDir + triggerType + "/" + "Target Constrained Vx (mm)").fill(tarVert.getPosition().x());
- aida.histogram1D(plotDir + triggerType + "/" + "Target Constrained Vy (mm)").fill(tarVert.getPosition().y());
- aida.histogram1D(plotDir + triggerType + "/" + "Target Constrained Vz (mm)").fill(tarVert.getPosition().z());
- aida.histogram1D(plotDir + triggerType + "/" + "Target Constrained Mass (GeV)").fill(tarV0.getMass());
- aida.histogram1D(plotDir + triggerType + "/" + "Target Constrained Chi2").fill(tarVert.getChi2());
+ tarconVx.fill(tarVert.getPosition().x());
+ tarconVy.fill(tarVert.getPosition().y());
+ tarconVz.fill(tarVert.getPosition().z());
+ tarconMass.fill(tarV0.getMass());
+ tarconChi2.fill(tarVert.getChi2());
}
List<ReconstructedParticle> finalStateParticles = event.get(ReconstructedParticle.class, finalStateParticlesColName);
- if (debug)
+ if (debug) {
System.out.println("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 (debug)
+ if (debug) {
System.out.println("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)
+ if (charge < 1) {
+ if (ele1 == null) {
ele1 = fsPart;
- else
+ } else if (!hasSharedStrips(ele1, fsPart, hitToStrips, hitToRotated)) {
ele2 = fsPart;
+ }
+ }
}
}
sumChargeHisto.fill(sumCharge);
@@ -339,8 +375,9 @@
@Override
public void printDQMData() {
System.out.println("V0Monitoring::printDQMData");
- for (Entry<String, Double> entry : monitoredQuantityMap.entrySet())
+ for (Entry<String, Double> entry : monitoredQuantityMap.entrySet()) {
System.out.println(entry.getKey() + " = " + entry.getValue());
+ }
System.out.println("*******************************");
}
@@ -353,9 +390,6 @@
IAnalysisFactory analysisFactory = IAnalysisFactory.create();
IFitFactory fitFactory = analysisFactory.createFitFactory();
IFitter fitter = fitFactory.createFitter("chi2");
- IHistogram1D bsconVx = aida.histogram1D(plotDir + triggerType + "/" + "BS Constrained Vx (mm)");
- IHistogram1D bsconVy = aida.histogram1D(plotDir + triggerType + "/" + "BS Constrained Vy (mm)");
- IHistogram1D bsconVz = aida.histogram1D(plotDir + triggerType + "/" + "BS Constrained Vz (mm)");
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};
@@ -368,8 +402,9 @@
double[] parsVy = resVy.fittedParameters();
double[] parsVz = resVz.fittedParameters();
- for (int i = 0; i < 5; i++)
+ for (int i = 0; i < 5; i++) {
System.out.println("Vertex Fit Parameters: " + resVx.fittedParameterNames()[i] + " = " + parsVx[i] + "; " + parsVy[i] + "; " + parsVz[i]);
+ }
IPlotter plotter = analysisFactory.createPlotterFactory().create("Vertex Position");
plotter.createRegions(1, 3);
@@ -383,12 +418,13 @@
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 (IOException ex) {
Logger.getLogger(V0Monitoring.class.getName()).log(Level.SEVERE, null, ex);
}
+ }
// monitoredQuantityMap.put(fpQuantNames[2], sumVx / nTotV0);
// monitoredQuantityMap.put(fpQuantNames[3], sumVy / nTotV0);
@@ -409,7 +445,9 @@
@Override
public void printDQMStrings() {
for (int i = 0; i < 9; i++)//TODO: do this in a smarter way...loop over the map
+ {
System.out.println("ALTER TABLE dqm ADD " + fpQuantNames[i] + " double;");
+ }
}
IFitResult fitVertexPosition(IHistogram1D h1d, IFitter fitter, double[] init, String range
@@ -423,14 +461,6 @@
return ifr;
}
- private double getMomentum(Track trk) {
-
- double px = trk.getTrackStates().get(0).getMomentum()[0];
- double py = trk.getTrackStates().get(0).getMomentum()[1];
- double pz = trk.getTrackStates().get(0).getMomentum()[2];
- return Math.sqrt(px * px + py * py + pz * pz);
- }
-
private BilliorVertex fitVertex(BilliorTrack electron, BilliorTrack positron) {
// Create a vertex fitter from the magnetic field.
double bField = 0.24;
@@ -454,4 +484,11 @@
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);
+ }
}
|
