Author: [log in to unmask] Date: Wed Jun 24 19:02:53 2015 New Revision: 3200 Log: add more DQM plots Modified: java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/TrackingMonitoring.java java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/TrackingResiduals.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 java/trunk/users/src/main/java/org/hps/users/meeg/ReconParticleCleanupDriver.java Modified: java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/TrackingMonitoring.java ============================================================================= --- java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/TrackingMonitoring.java (original) +++ java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/TrackingMonitoring.java Wed Jun 24 19:02:53 2015 @@ -73,6 +73,21 @@ IHistogram1D[] deltaTOnTrackTop = new IHistogram1D[nmodules]; IHistogram1D[] deltaTOnTrackBot = new IHistogram1D[nmodules]; + IHistogram1D trkChi2 = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Track Chi2", 50, 0, 25.0); + IHistogram1D nTracks = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Tracks per Event", 6, 0, 6); + IHistogram1D trkd0; + IHistogram1D trkphi; + IHistogram1D trkomega; + IHistogram1D trklam; + IHistogram1D trkz0; + IHistogram1D nHits; + IHistogram1D trackMeanTime; + IHistogram1D trackRMSTime; + IHistogram2D trackChi2RMSTime; + IHistogram1D seedRMSTime; + IHistogram2D trigTimeTrackTime; + IHistogram1D trigTime; + IHistogram1D trkYAtECALTop; IHistogram1D trkYAtECALBot; @@ -137,7 +152,7 @@ double omegaCut = 0.0005; double lambdaCut = 0.1; double z0Cut = 1.0; - double timeCut=24.0; + double timeCut = 24.0; public void setHelicalTrackHitCollectionName(String helicalTrackHitCollectionName) { this.helicalTrackHitCollectionName = helicalTrackHitCollectionName; @@ -152,18 +167,20 @@ this.detector = detector; aida.tree().cd("/"); - IHistogram1D trkChi2 = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Track Chi2", 50, 0, 25.0); - IHistogram1D nTracks = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Tracks per Event", 6, 0, 6); - IHistogram1D trkd0 = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "d0 ", 50, -d0Cut, d0Cut); - IHistogram1D trkphi = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "sinphi ", 50, -phiCut, phiCut); - IHistogram1D trkomega = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "omega ", 50, -omegaCut, omegaCut); - IHistogram1D trklam = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "tan(lambda) ", 50, -lambdaCut, lambdaCut); - IHistogram1D trkz0 = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "z0 ", 50, -z0Cut, z0Cut); - IHistogram1D nHits = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Hits per Track", 4, 3, 7); - IHistogram1D trackMeanTime = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Mean time of hits on track", 100, -timeCut, timeCut); - IHistogram1D trackRMSTime = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "RMS time of hits on track", 100, 0., 15.); - IHistogram2D trackChi2RMSTime = aida.histogram2D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Track chi2 vs. RMS time of hits", 100, 0., 15., 25, 0, 25.0); - IHistogram1D seedRMSTime = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "RMS time of hits on seed layers", 100, 0., 15.); + trkChi2 = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Track Chi2", 50, 0, 25.0); + nTracks = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Tracks per Event", 6, 0, 6); + trkd0 = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "d0 ", 50, -d0Cut, d0Cut); + trkphi = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "sinphi ", 50, -phiCut, phiCut); + trkomega = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "omega ", 50, -omegaCut, omegaCut); + trklam = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "tan(lambda) ", 50, -lambdaCut, lambdaCut); + trkz0 = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "z0 ", 50, -z0Cut, z0Cut); + nHits = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Hits per Track", 4, 3, 7); + trackMeanTime = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Mean time of hits on track", 100, -timeCut, timeCut); + trackRMSTime = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "RMS time of hits on track", 100, 0., 15.); + trackChi2RMSTime = aida.histogram2D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Track chi2 vs. RMS time of hits", 100, 0., 15., 25, 0, 25.0); + trigTimeTrackTime = aida.histogram2D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Trigger phase vs. mean time of hits", 100, -timeCut, timeCut, 6, 0, 24.0); + trigTime = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Trigger phase", 6, 0., 24.); + seedRMSTime = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "RMS time of hits on seed layers", 100, 0., 15.); trkYAtECALTop = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Track Y at ECAL: Top", 100, 0, 100); trkYAtECALBot = aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Track Y at ECAL: Bot", 100, 0, 100); for (int i = 1; i <= nmodules; i++) { @@ -253,50 +270,56 @@ aida.tree().cd("/"); - if (!event.hasCollection(LCRelation.class, helicalTrackHitRelationsCollectionName) || !event.hasCollection(LCRelation.class, rotatedHelicalTrackHitRelationsCollectionName)) + if (!event.hasCollection(LCRelation.class, helicalTrackHitRelationsCollectionName) || !event.hasCollection(LCRelation.class, rotatedHelicalTrackHitRelationsCollectionName)) { return; + } 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) + 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) + for (LCRelation relation : rotaterelations) { + if (relation != null && relation.getFrom() != null && relation.getTo() != null) { hittorotated.add(relation.getFrom(), relation.getTo()); - - if (!event.hasCollection(TrackerHit.class, helicalTrackHitCollectionName)) + } + } + + if (!event.hasCollection(TrackerHit.class, helicalTrackHitCollectionName)) { return; + } //check to see if this event is from the correct trigger (or "all"); - if (!matchTrigger(event)) + if (!matchTrigger(event)) { return; + } /* This doesn't work on reco'ed files...fix me! */ int[] topHits = {0, 0, 0, 0, 0, 0}; int[] botHits = {0, 0, 0, 0, 0, 0}; List<TrackerHit> hth = event.get(TrackerHit.class, helicalTrackHitCollectionName); for (TrackerHit hit : hth) { - int module = -99; int layer = ((RawTrackerHit) hit.getRawHits().get(0)).getLayerNumber(); - module = layer/2 + 1; + int module = layer / 2 + 1; Collection<TrackerHit> htsList = hittostrip.allFrom(hit); double hitTimes[] = new double[2]; for (TrackerHit hts : htsList) { int stripLayer = ((HpsSiSensor) ((RawTrackerHit) hts.getRawHits().get(0)).getDetectorElement()).getLayerNumber(); - hitTimes[stripLayer%2] = hts.getTime(); + hitTimes[stripLayer % 2] = hts.getTime(); } if (hit.getPosition()[1] > 0) { topHits[module - 1]++; xvsyTop[module - 1].fill(hit.getPosition()[0], hit.getPosition()[1]); - timevstimeTop[module - 1].fill(hitTimes[0],hitTimes[1]); + timevstimeTop[module - 1].fill(hitTimes[0], hitTimes[1]); } else { botHits[module - 1]++; xvsyBot[module - 1].fill(hit.getPosition()[0], Math.abs(hit.getPosition()[1])); - timevstimeBot[module - 1].fill(hitTimes[0],hitTimes[1]); + timevstimeBot[module - 1].fill(hitTimes[0], hitTimes[1]); } } @@ -306,13 +329,13 @@ } if (!event.hasCollection(Track.class, trackCollectionName)) { - aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Tracks per Event").fill(0); + nTracks.fill(0); return; } nEvents++; List<Track> tracks = event.get(Track.class, trackCollectionName); nTotTracks += tracks.size(); - aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Tracks per Event").fill(tracks.size()); + nTracks.fill(tracks.size()); int cntEle = 0; int cntPos = 0; int cntTop = 0; @@ -321,10 +344,11 @@ for (Track trk : tracks) { Hep3Vector trackPosAtEcalFace = TrackUtils.extrapolateTrack(trk, ecalFace); double yAtECal = trackPosAtEcalFace.y(); - if (yAtECal > 0) + if (yAtECal > 0) { trkYAtECALTop.fill(yAtECal); - else + } else { trkYAtECALBot.fill(Math.abs(yAtECal)); + } nTotHits += trk.getTrackerHits().size(); double d0 = trk.getTrackStates().get(0).getD0(); @@ -332,13 +356,13 @@ double omega = trk.getTrackStates().get(0).getOmega(); double lambda = trk.getTrackStates().get(0).getTanLambda(); double z0 = trk.getTrackStates().get(0).getZ0(); - aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Track Chi2").fill(trk.getChi2()); - aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Hits per Track").fill(trk.getTrackerHits().size()); - aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "d0 ").fill(trk.getTrackStates().get(0).getD0()); - aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "sinphi ").fill(Math.sin(trk.getTrackStates().get(0).getPhi())); - aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "omega ").fill(trk.getTrackStates().get(0).getOmega()); - aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "tan(lambda) ").fill(trk.getTrackStates().get(0).getTanLambda()); - aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "z0 ").fill(trk.getTrackStates().get(0).getZ0()); + trkChi2.fill(trk.getChi2()); + nHits.fill(trk.getTrackerHits().size()); + trkd0.fill(trk.getTrackStates().get(0).getD0()); + trkphi.fill(Math.sin(trk.getTrackStates().get(0).getPhi())); + trkomega.fill(trk.getTrackStates().get(0).getOmega()); + trklam.fill(trk.getTrackStates().get(0).getTanLambda()); + trkz0.fill(trk.getTrackStates().get(0).getZ0()); d0VsPhi0.fill(d0, sinphi0); d0Vsomega.fill(d0, omega); d0Vslambda.fill(d0, lambda); @@ -393,19 +417,22 @@ rmsTime += Math.pow(hts.getTime() - meanTime, 2); HpsSiSensor sensor = (HpsSiSensor) ((RawTrackerHit) hts.getRawHits().get(0)).getDetectorElement(); int layer = sensor.getLayerNumber(); - if (layer <= 6) + if (layer <= 6) { rmsSeedTime += Math.pow(hts.getTime() - meanSeedTime, 2); + } String sensorName = getNiceSensorName(sensor); getSensorPlot(plotDir + trackCollectionName + "/" + triggerType + "/" + timeresidDir + "hitTimeResidual_", sensorName).fill((hts.getTime() - meanTime) * nStrips / (nStrips - 1)); //correct residual for bias } } rmsTime = Math.sqrt(rmsTime / nStrips); - aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Mean time of hits on track").fill(meanTime); - aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "RMS time of hits on track").fill(rmsTime); - aida.histogram2D(plotDir + trackCollectionName + "/" + triggerType + "/" + "Track chi2 vs. RMS time of hits").fill(rmsTime, trk.getChi2()); + trackMeanTime.fill(meanTime); + trackRMSTime.fill(rmsTime); + trackChi2RMSTime.fill(rmsTime, trk.getChi2()); + trigTimeTrackTime.fill(meanTime, event.getTimeStamp() % 24); + trigTime.fill(event.getTimeStamp() % 24); rmsSeedTime = Math.sqrt(rmsSeedTime / nSeedStrips); - aida.histogram1D(plotDir + trackCollectionName + "/" + triggerType + "/" + "RMS time of hits on seed layers").fill(rmsSeedTime); + seedRMSTime.fill(rmsSeedTime); if (trk.getTrackStates().get(0).getOmega() < 0) {//positrons trkChi2Pos.fill(trk.getChi2()); @@ -475,8 +502,9 @@ //IHistogram1D occupancyPlot = aida.histogram1D(sensor.getName().replaceAll("Tracker_TestRunModule_", ""), 640, 0, 639); IHistogram1D hitTimeResidual = getSensorPlot(plotDir + trackCollectionName + "/" + triggerType + "/" + timeresidDir + "hitTimeResidual_", getNiceSensorName(sensor)); IFitResult result = fitGaussian(hitTimeResidual, fitter, "range=\"(-20.0,20.0)\""); - if (result != null) + if (result != null) { System.out.format("%s\t%f\t%f\t%d\t%d\n", getNiceSensorName(sensor), result.fittedParameters()[1], result.fittedParameters()[2], sensor.getFebID(), sensor.getFebHybridID()); + } } monitoredQuantityMap.put(trackCollectionName + " " + triggerType + " " + trackingQuantNames[0], (double) nTotTracks / nEvents); @@ -503,15 +531,17 @@ @Override public void printDQMData() { System.out.println("ReconMonitoring::printDQMData"); - for (Map.Entry<String, Double> entry : monitoredQuantityMap.entrySet()) + for (Map.Entry<String, Double> entry : monitoredQuantityMap.entrySet()) { System.out.println(entry.getKey() + " = " + entry.getValue()); + } System.out.println("*******************************"); } @Override public void printDQMStrings() { - for (Map.Entry<String, Double> entry : monitoredQuantityMap.entrySet()) + for (Map.Entry<String, Double> entry : monitoredQuantityMap.entrySet()) { System.out.println("ALTER TABLE dqm ADD " + entry.getKey() + " double;"); + } } private IHistogram1D getSensorPlot(String prefix, HpsSiSensor sensor) { Modified: java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/TrackingResiduals.java ============================================================================= --- java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/TrackingResiduals.java (original) +++ java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/TrackingResiduals.java Wed Jun 24 19:02:53 2015 @@ -71,8 +71,6 @@ @Override public void process(EventHeader event) { aida.tree().cd("/"); - if (!event.hasCollection(GenericObject.class, trackTimeDataCollectionName)) - return; if (!event.hasCollection(GenericObject.class, trackResidualsCollectionName)) return; //check to see if this event is from the correct trigger (or "all"); @@ -94,12 +92,16 @@ } } - List<GenericObject> ttdList = event.get(GenericObject.class, trackTimeDataCollectionName); - for (GenericObject ttd : ttdList) { - int nResid = ttd.getNDouble(); - for (int i = 1; i <= nResid; i++) - aida.histogram1D(plotDir + triggerType + "/"+timeresDir + "HalfModule " + i + " t Residual").fill(ttd.getDoubleVal(i - 1));//x is the double value in the generic object - } + if (event.hasCollection(GenericObject.class, trackTimeDataCollectionName)) { + List<GenericObject> ttdList = event.get(GenericObject.class, trackTimeDataCollectionName); + for (GenericObject ttd : ttdList) { + int nResid = ttd.getNDouble(); + for (int i = 1; i <= nResid; i++) { + aida.histogram1D(plotDir + triggerType + "/" + timeresDir + "HalfModule " + i + " t Residual").fill(ttd.getDoubleVal(i - 1));//x is the double value in the generic object + } + } + } + if (!event.hasCollection(GenericObject.class, gblStripClusterDataCollectionName)) return; List<GenericObject> gblSCDList = event.get(GenericObject.class, gblStripClusterDataCollectionName); 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 Wed Jun 24 19:02:53 2015 @@ -48,6 +48,14 @@ IHistogram1D vertexX; IHistogram1D vertexY; IHistogram1D vertexZ; + IHistogram1D vertexPx; + IHistogram1D vertexPy; + IHistogram1D vertexPz; + IHistogram1D vertexU; + IHistogram1D vertexV; +// IHistogram1D vertexW; + IHistogram2D vertexVZ; + IHistogram2D vertexZY; IHistogram1D deltaP; IHistogram1D deltaPRad; @@ -129,6 +137,14 @@ vertexX = aida.histogram1D(plotDir + triggerType + "/" + "Vertex X Position (mm)", 100, -v0VxMax, v0VxMax); vertexY = aida.histogram1D(plotDir + triggerType + "/" + "Vertex Y Position (mm)", 100, -v0VyMax, v0VyMax); vertexZ = aida.histogram1D(plotDir + triggerType + "/" + "Vertex Z Position (mm)", 100, -v0VzMax, v0VzMax); + vertexPx = aida.histogram1D(plotDir + triggerType + "/" + "Vertex Px (GeV)", 100, -0.1, 0.1); + vertexPy = aida.histogram1D(plotDir + triggerType + "/" + "Vertex Py (GeV)", 100, -0.1, 0.1); + vertexPz = aida.histogram1D(plotDir + triggerType + "/" + "Vertex Pz (GeV)", 100, 0.0, v0PzMax); + vertexU = aida.histogram1D(plotDir + triggerType + "/" + "Vertex Px over Ptot (GeV)", 100, -0.1, 0.1); + vertexV = aida.histogram1D(plotDir + triggerType + "/" + "Vertex Py over Ptot (GeV)", 100, -0.1, 0.1); +// vertexW = aida.histogram1D(plotDir + triggerType + "/" + "Vertex Pz overPtot (GeV)", 100, 0.95, 1.0); + vertexVZ = aida.histogram2D(plotDir + triggerType + "/" + "Vertex Py over Ptot vs. Z", 100, -v0VzMax, v0VzMax, 100, -0.1, 0.1); + vertexZY = aida.histogram2D(plotDir + triggerType + "/" + "Vertex Z vs. Y", 100, -v0VyMax, v0VyMax, 100, -v0VzMax, v0VzMax); } @Override @@ -185,23 +201,23 @@ // 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; + continue; } if (Math.abs(v0Mom.y()) > v0PyMax) { - break; + continue; } if (Math.abs(v0Mom.x()) > v0PxMax) { - break; + continue; } Hep3Vector v0Vtx = uncVert.getPosition(); if (Math.abs(v0Vtx.z()) > v0VzMax) { - break; + continue; } if (Math.abs(v0Vtx.y()) > v0VyMax) { - break; + continue; } if (Math.abs(v0Vtx.x()) > v0VxMax) { - break; + continue; } List<Track> tracks = new ArrayList<Track>(); @@ -233,6 +249,7 @@ boolean trackTimeDiffCut = Math.abs(trackTimes.get(0) - trackTimes.get(1)) < trkTimeDiff; boolean pCut = electron.getMomentum().magnitude() > trkPzMin && positron.getMomentum().magnitude() > trkPzMin; boolean pTotCut = uncV0.getMomentum().magnitude() > v0PzMin && uncV0.getMomentum().magnitude() < v0PzMax; +// if (electron.getMomentum().y()* positron.getMomentum().y()<0) continue; if (trackTimeDiffCut) { vertexMassMomentum.fill(uncV0.getMomentum().magnitude(), uncV0.getMass()); vertexedTrackMomentum2D.fill(electron.getMomentum().magnitude(), positron.getMomentum().magnitude()); @@ -250,6 +267,14 @@ vertexX.fill(v0Vtx.x()); vertexY.fill(v0Vtx.y()); vertexZ.fill(v0Vtx.z()); + vertexPx.fill(uncV0.getMomentum().x()); + vertexPy.fill(uncV0.getMomentum().y()); + vertexPz.fill(uncV0.getMomentum().z()); + vertexU.fill(uncV0.getMomentum().x()/uncV0.getMomentum().magnitude()); + vertexV.fill(uncV0.getMomentum().y()/uncV0.getMomentum().magnitude()); +// vertexW.fill(uncV0.getMomentum().z()/uncV0.getMomentum().magnitude()); + vertexVZ.fill(v0Vtx.z(), uncV0.getMomentum().y()/uncV0.getMomentum().magnitude()); + vertexZY.fill(v0Vtx.y(), v0Vtx.z()); } } // System.out.println(tarV0.getTracks()) 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 Wed Jun 24 19:02:53 2015 @@ -65,6 +65,12 @@ /* beamspot constrained */ IHistogram1D nV0; + + IHistogram1D v0Time; + IHistogram1D v0Dt; + IHistogram2D trigTimeV0Time; + IHistogram1D trigTime; + IHistogram1D bsconMass; IHistogram1D bsconVx; IHistogram1D bsconVy; @@ -144,6 +150,11 @@ /* beamspot constrained */ nV0 = aida.histogram1D(plotDir + triggerType + "/" + "Number of V0 per event", 10, 0, 10); + v0Time = aida.histogram1D(plotDir + triggerType + "/" + "V0 mean time", 100, -25, 25); + v0Dt = aida.histogram1D(plotDir + triggerType + "/" + "V0 time difference", 100, -25, 25); + trigTimeV0Time = aida.histogram2D(plotDir + triggerType + "/" + "Trigger phase vs. V0 mean time", 100, -25, 25, 6, 0, 24); + trigTime = aida.histogram1D(plotDir + triggerType + "/" + "Trigger phase", 6, 0, 24); + 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); @@ -259,6 +270,14 @@ { 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); Modified: java/trunk/users/src/main/java/org/hps/users/meeg/ReconParticleCleanupDriver.java ============================================================================= --- java/trunk/users/src/main/java/org/hps/users/meeg/ReconParticleCleanupDriver.java (original) +++ java/trunk/users/src/main/java/org/hps/users/meeg/ReconParticleCleanupDriver.java Wed Jun 24 19:02:53 2015 @@ -12,7 +12,7 @@ /** * Remove final state particles with bad track-cluster time matching, and - * vertices with shared tracks. + * vertices with shared hits. * * @author Sho Uemura <[log in to unmask]> * @version $Id: $ @@ -20,11 +20,16 @@ public class ReconParticleCleanupDriver extends Driver { private final String finalStateParticlesColName = "FinalStateParticles"; - private final String unconstrainedV0CandidatesColName = "UnconstrainedV0Candidates"; - private final String beamConV0CandidatesColName = "BeamspotConstrainedV0Candidates"; - private final String targetV0ConCandidatesColName = "TargetConstrainedV0Candidates"; + private final String[] v0ColNames = {"UnconstrainedV0Candidates", "BeamspotConstrainedV0Candidates", "TargetConstrainedV0Candidates"}; + private double fsDeltaT = 43.5; - private double fsDeltaTCut = 5.0; + private double fsDeltaTCut = -1; + private double maxTrackDt = -1; + private boolean discardUnmatchedTracks = false; + + public void setDiscardUnmatchedTracks(boolean discardUnmatchedTracks) { + this.discardUnmatchedTracks = discardUnmatchedTracks; + } /** * Center value for (cluster - track) time cut. @@ -36,12 +41,23 @@ } /** - * Cut window half-width for (cluster - track) time cut. + * Cut window half-width for (cluster - track) time cut. Negative value + * disables this cut. * * @param fsDeltaTCut */ public void setFsDeltaTCut(double fsDeltaTCut) { this.fsDeltaTCut = fsDeltaTCut; + } + + /** + * Cut window half-width for (track - track) time cut. Negative value + * disables this cut. + * + * @param maxTrackDt + */ + public void setMaxTrackDt(double maxTrackDt) { + this.maxTrackDt = maxTrackDt; } @Override @@ -51,7 +67,7 @@ for (Iterator<ReconstructedParticle> iter = event.get(ReconstructedParticle.class, finalStateParticlesColName).listIterator(); iter.hasNext();) { ReconstructedParticle fs = iter.next(); - if (fs.getClusters().isEmpty()) {//track without cluster, discard + if (discardUnmatchedTracks && fs.getClusters().isEmpty()) {//track without cluster, discard iter.remove(); continue; } @@ -60,55 +76,50 @@ continue; } - double deltaT = ClusterUtilities.getSeedHitTime(fs.getClusters().get(0)) - TrackUtils.getTrackTime(fs.getTracks().get(0), hitToStrips, hitToRotated); - if (Math.abs(deltaT - fsDeltaT) > fsDeltaTCut) {//bad track-cluster time match, discard - iter.remove(); + if (!fs.getClusters().isEmpty() && !fs.getTracks().isEmpty()) { + double deltaT = ClusterUtilities.getSeedHitTime(fs.getClusters().get(0)) - TrackUtils.getTrackTime(fs.getTracks().get(0), hitToStrips, hitToRotated); + if (fsDeltaTCut > 0 && Math.abs(deltaT - fsDeltaT) > fsDeltaTCut) {//bad track-cluster time match, discard + iter.remove(); + } } } Set<ReconstructedParticle> fsParticles = new HashSet<ReconstructedParticle>(event.get(ReconstructedParticle.class, finalStateParticlesColName)); - v0Loop: - for (Iterator<ReconstructedParticle> iter = event.get(ReconstructedParticle.class, unconstrainedV0CandidatesColName).listIterator(); iter.hasNext();) { - ReconstructedParticle v0 = iter.next(); - if (hasSharedStrips(v0, hitToStrips, hitToRotated)) { - iter.remove(); - continue; - } - for (ReconstructedParticle particle : v0.getParticles()) { - if (!fsParticles.contains(particle)) { + for (String colName : v0ColNames) { + v0Loop: + for (Iterator<ReconstructedParticle> iter = event.get(ReconstructedParticle.class, colName).listIterator(); iter.hasNext();) { + ReconstructedParticle v0 = iter.next(); + + ReconstructedParticle[] particles = new ReconstructedParticle[2]; + for (ReconstructedParticle particle : v0.getParticles()) // tracks.addAll(particle.getTracks()); //add add electron first, then positron...down below + { + if (particle.getCharge() < 0) { + particles[0] = particle; + } else if (particle.getCharge() > 0) { + particles[1] = particle; + } else { + throw new RuntimeException("expected only electron and positron in vertex, got something with charge 0"); + } + } + if (particles[0] == null || particles[1] == null) { + throw new RuntimeException("vertex needs e+ and e- but is missing one or both"); + } + double deltaT = TrackUtils.getTrackTime(particles[0].getTracks().get(0), hitToStrips, hitToRotated) - TrackUtils.getTrackTime(particles[1].getTracks().get(0), hitToStrips, hitToRotated); //electron time - positron time + + if (hasSharedStrips(v0, hitToStrips, hitToRotated)) { iter.remove(); - continue v0Loop; + continue; } - } - } - - v0Loop: - for (Iterator<ReconstructedParticle> iter = event.get(ReconstructedParticle.class, beamConV0CandidatesColName).listIterator(); iter.hasNext();) { - ReconstructedParticle v0 = iter.next(); - if (hasSharedStrips(v0, hitToStrips, hitToRotated)) { - iter.remove(); - continue; - } - for (ReconstructedParticle particle : v0.getParticles()) { - if (!fsParticles.contains(particle)) { + if (maxTrackDt > 0 && Math.abs(deltaT) > maxTrackDt) { iter.remove(); - continue v0Loop; + continue; } - } - } - - v0Loop: - for (Iterator<ReconstructedParticle> iter = event.get(ReconstructedParticle.class, targetV0ConCandidatesColName).listIterator(); iter.hasNext();) { - ReconstructedParticle v0 = iter.next(); - if (hasSharedStrips(v0, hitToStrips, hitToRotated)) { - iter.remove(); - continue; - } - for (ReconstructedParticle particle : v0.getParticles()) { - if (!fsParticles.contains(particle)) { - iter.remove(); - continue v0Loop; + for (ReconstructedParticle particle : v0.getParticles()) { + if (!fsParticles.contains(particle)) { + iter.remove(); + continue v0Loop; + } } } } @@ -121,5 +132,4 @@ 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); } - }