Commit in hps-java/src/main/java/org/lcsim/hps on MAIN | |||
users/mgraham/TestAnalysisDriver.java | +114 | -36 | 1.1 -> 1.2 |
recon/tracking/DumbShaperFit.java | +11 | -2 | 1.1 -> 1.2 |
/HPSStripMaker.java | +1 | -1 | 1.3 -> 1.4 |
/HPSShapeFitParameters.java | +12 | -2 | 1.3 -> 1.4 |
/HPSNearestNeighborRMS.java | +14 | -3 | 1.3 -> 1.4 |
/HPSRawTrackerHitFitterDriver.java | +1 | -1 | 1.5 -> 1.6 |
/HPSSVTCalibrationConstants.java | +2 | -2 | 1.3 -> 1.4 |
/HPSShaperAnalyticFitAlgorithm.java | +8 | -8 | 1.3 -> 1.4 |
+163 | -55 |
Very minor tweaks.
diff -u -r1.1 -r1.2 --- TestAnalysisDriver.java 25 Apr 2012 15:23:42 -0000 1.1 +++ TestAnalysisDriver.java 25 Apr 2012 18:01:32 -0000 1.2 @@ -4,53 +4,131 @@
*/ package org.lcsim.hps.users.mgraham;
+import hep.aida.IAnalysisFactory; +import hep.aida.IProfile1D; +import java.io.IOException; +import java.util.HashMap;
import java.util.List;
-import org.lcsim.event.EventHeader; -import org.lcsim.event.LCRelation; -import org.lcsim.event.RelationalTable; -import org.lcsim.event.Track;
+import java.util.Map; +import java.util.logging.Level; +import java.util.logging.Logger; +import org.lcsim.event.*;
import org.lcsim.event.base.BaseRelationalTable;
+import org.lcsim.fit.helicaltrack.HelicalTrackFit; +import org.lcsim.fit.helicaltrack.HelicalTrackHit; +import org.lcsim.hps.recon.tracking.HPSFittedRawTrackerHit;
import org.lcsim.hps.recon.tracking.TrackAnalysis;
+import org.lcsim.recon.tracking.digitization.sisim.SiTrackerHitStrip1D; +import org.lcsim.recon.tracking.seedtracker.SeedCandidate; +import org.lcsim.recon.tracking.seedtracker.SeedTrack;
import org.lcsim.util.Driver;
+import org.lcsim.util.aida.AIDA;
/** * * @author mgraham */
-public class TestAnalysisDriver extends Driver{ - int nevents=0; - int naccepted=0; - - public void process(
+public class TestAnalysisDriver extends Driver { + + int nevents = 0; + int naccepted = 0; + private AIDA aida = AIDA.defaultInstance(); + private IAnalysisFactory af = aida.analysisFactory(); + public String outputPlots = "myplots.aida"; + Map<String, IProfile1D> clsizeMap = new HashMap<String, IProfile1D>(); + String[] detNames = {"Tracker"}; + Integer[] nlayers = {10}; + + public void process(
EventHeader event) {
- nevents++; - List<Track> tracklist = event.get(Track.class, "MatchedTracks"); - if(tracklist.size()<2) - return; - - RelationalTable hittomc = new BaseRelationalTable(RelationalTable.Mode.MANY_TO_MANY, RelationalTable.Weighting.UNWEIGHTED); - List<LCRelation> mcrelations = event.get(LCRelation.class, "HelicalTrackMCRelations");
+ nevents++;
- for (LCRelation relation : mcrelations) { - if (relation != null && relation.getFrom() != null && relation.getTo() != null) { - hittomc.add(relation.getFrom(), relation.getTo()); - } - } - int ok=0; - for (Track track : tracklist) { //remember, these tracks are in the lcsim tracking frame! - TrackAnalysis tkanal = new TrackAnalysis(track, hittomc); - if(Math.abs(tkanal.getMCParticle().getPDGID())==611) - ok++; - //do some stuff to makes sure tracks are great - //is there an e+e- from the muonium
+ List<RawTrackerHit> rawHits = event.get(RawTrackerHit.class, "SVTRawTrackerHits"); + List<HPSFittedRawTrackerHit> fittedrawHits = event.get(HPSFittedRawTrackerHit.class, "SVTFittedRawTrackerHits"); + + List<SiTrackerHitStrip1D> stripHits = event.get(SiTrackerHitStrip1D.class, "StripClusterer_SiTrackerHitStrip1D"); + + + +Map<Track, TrackAnalysis> tkanalMap = new HashMap<Track, TrackAnalysis>(); + RelationalTable nearestHitToTrack = new BaseRelationalTable(RelationalTable.Mode.MANY_TO_MANY, RelationalTable.Weighting.UNWEIGHTED); + Map<Track, Double> l1Isolation = new HashMap<Track, Double>(); + Map<Track, Double> l1DeltaZ = new HashMap<Track, Double>(); + + String fittedHitsDir = "FittedHits/"; + for(HPSFittedRawTrackerHit hrth:fittedrawHits){ + double fittedAmp=hrth.getAmp(); + double fittedT0=hrth.getT0(); + String sensorName=hrth.getRawTrackerHit().getDetectorElement().getName(); + aida.cloud1D(fittedHitsDir + sensorName+" fitted Amplitude").fill(fittedAmp); + aida.cloud1D(fittedHitsDir + sensorName+"fitted T0").fill(fittedT0);
}
- if(ok==2) - naccepted++; - - } - public void endOfData() { - - System.out.println("# of muonium events= " + naccepted + "; # of total = " + nevents);
+ + /* + List<Track> tracklist = event.get(Track.class, "MatchedTracks"); + List<HelicalTrackHit> hthits = event.get(HelicalTrackHit.class, "HelicalTrackHits"); + aida.cloud1D("HelicalTrackHits per Event").fill(hthits.size()); + aida.cloud1D("Matched Tracks per Event").fill(tracklist.size()); + String trackdir = "TrackInfo/"; + // Analyze the tracks in the event + for (Track track : tracklist) { + // Calculate the track pT and cos(theta) + double d0 = track.getTrackParameter(HelicalTrackFit.dcaIndex); + double z0 = track.getTrackParameter(HelicalTrackFit.z0Index); + double phi0 = track.getTrackParameter(HelicalTrackFit.phi0Index); + double slope = track.getTrackParameter(HelicalTrackFit.slopeIndex); + double curve = track.getTrackParameter(HelicalTrackFit.curvatureIndex); + double d0Err = Math.sqrt(track.getErrorMatrix().e(HelicalTrackFit.dcaIndex, HelicalTrackFit.dcaIndex)); + double z0Err = Math.sqrt(track.getErrorMatrix().e(HelicalTrackFit.z0Index, HelicalTrackFit.z0Index)); + double phi0Err = Math.sqrt(track.getErrorMatrix().e(HelicalTrackFit.phi0Index, HelicalTrackFit.phi0Index)); + double slopeErr = Math.sqrt(track.getErrorMatrix().e(HelicalTrackFit.slopeIndex, HelicalTrackFit.slopeIndex)); + double curveErr = Math.sqrt(track.getErrorMatrix().e(HelicalTrackFit.curvatureIndex, HelicalTrackFit.curvatureIndex)); + SeedTrack stEle = (SeedTrack) track; + SeedCandidate seedEle = stEle.getSeedCandidate(); + HelicalTrackFit ht = seedEle.getHelix(); + + double xoca = ht.x0(); + double yoca = ht.y0(); + double[] poca = {xoca, yoca, z0}; + double mom[] = track.getMomentum(); + double px = mom[0]; + double py = mom[1]; + double pz = mom[2]; + double pperp = Math.sqrt(py * py + pz * pz); + double pt = Math.sqrt(px * px + py * py); + double p = Math.sqrt(pt * pt + pz * pz); + double phi = Math.atan2(py, px); + double cth = pz / Math.sqrt(pt * pt + pz * pz); + double sth = pt / Math.sqrt(pt * pt + pz * pz); + double th = Math.atan2(pt, pz); + fillTrackInfo(trackdir, "all tracks", track.getChi2(), p, pperp, px, py, pz, phi, cth, d0, xoca, yoca, z0); + } + */ + } + private void fillTrackInfo(String dir, String species, double chi2, double p, double pperp, double px, double py, double pz, double phi, double cth, double doca, double xoca, double yoca, double zoca) { + aida.cloud1D(dir + "total chi^2 for " + species).fill(chi2); + +// aida.cloud1D(trackdir + "circle chi^2 for " + species).fill(ht.chisq()[0]); +// aida.cloud1D(trackdir + "linear chi^2 for " + species).fill(ht.chisq()[1] + aida.cloud1D(dir + "p for " + species).fill(p); + aida.cloud1D(dir + "pperp for " + species).fill(pperp); + aida.cloud1D(dir + "px for " + species).fill(px); + aida.cloud1D(dir + "py for " + species).fill(py); + aida.cloud1D(dir + "pz for " + species).fill(pz); + aida.cloud1D(dir + "phi for " + species).fill(phi); + aida.cloud1D(dir + "cos(theta) for " + species).fill(cth); + aida.cloud1D(dir + "DOCA for " + species).fill(doca); + aida.cloud1D(dir + "XOCA for " + species).fill(xoca); + aida.cloud1D(dir + "YOCA for " + species).fill(yoca); + aida.cloud1D(dir + "ZOCA for " + species).fill(zoca); + aida.cloud2D(dir + "doca vs xoca for " + species).fill(xoca, doca); + } + public void endOfData() { + try { + aida.saveAs(outputPlots); + } catch (IOException ex) { + Logger.getLogger(DetailedAnalysisDriver.class.getName()).log(Level.SEVERE, null, ex); + }
}
-}
+}
diff -u -r1.1 -r1.2 --- DumbShaperFit.java 25 Apr 2012 15:23:42 -0000 1.1 +++ DumbShaperFit.java 25 Apr 2012 18:01:32 -0000 1.2 @@ -15,8 +15,17 @@
@Override public HPSShapeFitParameters fitShape(RawTrackerHit rth, ChannelConstants constants) { HPSShapeFitParameters fitresults = new HPSShapeFitParameters();
- fitresults.setAmp(1000.666); - fitresults.setT0(5.999);
+ short[] adcVals=rth.getADCValues(); + double maxADC=-99999; + double t0=-999; + for(int i=0;i<6;i++){ + if(adcVals[i]>maxADC){ + maxADC=adcVals[i]; + t0=24.0*i; + } + } + fitresults.setAmp(maxADC); + fitresults.setT0(t0);
return fitresults; } }
diff -u -r1.3 -r1.4 --- HPSStripMaker.java 25 Apr 2012 15:23:42 -0000 1.3 +++ HPSStripMaker.java 25 Apr 2012 18:01:32 -0000 1.4 @@ -167,7 +167,7 @@
_max_cluster_nstrips = max_cluster_nstrips; }
- private SiTrackerHitStrip1D makeTrackerHit(List<HPSFittedRawTrackerHit> cluster, SiSensorElectrodes electrodes) {
+ private SiTrackerHitStrip1D makeTrackerHit(List<HPSFittedRawTrackerHit> cluster, SiSensorElectrodes electrodes) {
Hep3Vector position = getPosition(cluster, electrodes); SymmetricMatrix covariance = getCovariance(cluster, electrodes); double time = getTime(cluster);
diff -u -r1.3 -r1.4 --- HPSShapeFitParameters.java 24 Apr 2012 23:32:31 -0000 1.3 +++ HPSShapeFitParameters.java 25 Apr 2012 18:01:32 -0000 1.4 @@ -18,6 +18,7 @@
private double _ampErr = Double.NaN; private double _tp = Double.NaN; private double _tpErr = Double.NaN;
+ private double _chiSq=Double.NaN;
public HPSShapeFitParameters() { }
@@ -51,6 +52,10 @@
this._tpErr = _tpErr; }
+ public void setChiSq(double _chiSq) { + this._chiSq = _chiSq; + } +
public double getT0() { return _t0; }
@@ -74,6 +79,9 @@
public double getTpErr() { return _tpErr; }
+ public double getChiSq() { + return _chiSq; + }
@Override public int getNInt() {
@@ -87,7 +95,7 @@
@Override public int getNDouble() {
- return 6;
+ return 7;
} @Override
@@ -115,8 +123,10 @@
return _tp; case 5: return _tpErr;
+ case 6: + return _chiSq;
default:
- throw new UnsupportedOperationException("Only 6 double values in " + this.getClass().getSimpleName());
+ throw new UnsupportedOperationException("Only 7 double values in " + this.getClass().getSimpleName());
} }
diff -u -r1.3 -r1.4 --- HPSNearestNeighborRMS.java 24 Apr 2012 17:46:01 -0000 1.3 +++ HPSNearestNeighborRMS.java 25 Apr 2012 18:01:32 -0000 1.4 @@ -16,7 +16,9 @@
private static String _NAME = "NearestNeighborRMS"; private double _seed_threshold; private double _neighbor_threshold;
- private double _cluster_threshold;
+ private double _cluster_threshold; + private double _meanTime=48; + private double _timeWindow = 12;
/** * Instantiate NearestNeighborRMS with specified thresholds.
@@ -159,12 +161,12 @@
// Get the signal from the readout chip double signal = base_hit.getAmp(); double noiseRMS = HPSSVTCalibrationConstants.getNoise((SiSensor)rth.getDetectorElement(),channel_number);
-
+ double time=base_hit.getT0();
// Mark this hit as available for clustering if it is above the neighbor threshold clusterable.put(channel_number, signal/noiseRMS >= _neighbor_threshold); // Add this hit to the list of seeds if it is above the seed threshold
- if (signal/noiseRMS >= _seed_threshold) {
+ if (signal/noiseRMS >= _seed_threshold&&passTimingCut(base_hit)) {
cluster_seeds.add(channel_number); } }
@@ -236,6 +238,15 @@
return cluster_list; }
+ private boolean passTimingCut(HPSFittedRawTrackerHit hit){ + + boolean pass=false; + double time=hit.getT0(); + if((time-_meanTime)<_timeWindow) + pass = true; + + return pass; + }
public int getNeighborCell(int cell, int ncells_0, int ncells_1) {
diff -u -r1.5 -r1.6 --- HPSRawTrackerHitFitterDriver.java 25 Apr 2012 15:23:42 -0000 1.5 +++ HPSRawTrackerHitFitterDriver.java 25 Apr 2012 18:01:32 -0000 1.6 @@ -53,7 +53,7 @@
public void process(EventHeader event) { List<RawTrackerHit> rawHits = event.get(RawTrackerHit.class, rawHitCollectionName); if (rawHits == null) {
- throw new RuntimeException("Event is missing ECal hits collection!");
+ throw new RuntimeException("Event is missing Raw hits collection!");
} List<HPSFittedRawTrackerHit> hits = new ArrayList<HPSFittedRawTrackerHit>(); List<HPSShapeFitParameters> fits = new ArrayList<HPSShapeFitParameters>();
diff -u -r1.3 -r1.4 --- HPSSVTCalibrationConstants.java 25 Apr 2012 15:23:42 -0000 1.3 +++ HPSSVTCalibrationConstants.java 25 Apr 2012 18:01:32 -0000 1.4 @@ -105,8 +105,8 @@
Pair<SiSensor, Integer> sensorChannel = new Pair(sensor, channel); ChannelConstants constants = new ChannelConstants(); // ...don't have a constants file yet!
- noise.put(sensorChannel,40.0); - pedestal.put(sensorChannel,7000.0);
+ noise.put(sensorChannel,20.0); + pedestal.put(sensorChannel,1638.0);
tShaping.put(sensorChannel,36.0); //// constants.setNoise(noise.get(sensorChannel));
diff -u -r1.3 -r1.4 --- HPSShaperAnalyticFitAlgorithm.java 25 Apr 2012 16:23:03 -0000 1.3 +++ HPSShaperAnalyticFitAlgorithm.java 25 Apr 2012 18:01:32 -0000 1.4 @@ -6,7 +6,7 @@
/** * Fast fitter; currently only fits single hits. Uses Tp from ChannelConstants; fits values and errors for T0 and amplitude. * @author meeg
- * @version $Id: HPSShaperAnalyticFitAlgorithm.java,v 1.3 2012/04/25 16:23:03 meeg Exp $
+ * @version $Id: HPSShaperAnalyticFitAlgorithm.java,v 1.4 2012/04/25 18:01:32 mgraham Exp $
*/ public class HPSShaperAnalyticFitAlgorithm implements HPSShaperFitAlgorithm {
@@ -21,8 +21,8 @@
bestStart = i; } }
- fitSection(rth, constants, fit, bestStart); - return fit;
+ fit.setChiSq(fitSection(rth, constants, fit, bestStart)); + return fit;
} private double fitSection(RawTrackerHit rth, ChannelConstants constants, HPSShapeFitParameters fit, int start) {
@@ -30,9 +30,9 @@
double[] y = new double[length]; double[] t = new double[length];
- for (int i = 0; i < length; i++) { - y[i] = rth.getADCValues()[start + i] - constants.getPedestal(); - t[i] = HPSSVTConstants.SAMPLE_INTERVAL * i;
+ for (int i = 0; i < length; i++) { + y[i] = rth.getADCValues()[start + i] - constants.getPedestal(); + t[i] = HPSSVTConstants.SAMPLE_INTERVAL * i;
} double[] p = new double[length];
@@ -74,12 +74,12 @@
fit.setT0(t0); fit.setT0Err(Math.sqrt(time_var)); fit.setTp(constants.getTp());
-
+
double chisq = 0; for (int i = 0; i < length; i++) { double fit_y = A * (t[i] / constants.getTp()) * Math.exp(1 - t[i] / constants.getTp()); chisq += Math.pow((fit_y - y[i]) / constants.getNoise(), 2);
- }
+ }
return chisq / (length - 2); //TODO: p-value would be better here }
Use REPLY-ALL to reply to list
To unsubscribe from the LCD-CVS list, click the following link:
https://listserv.slac.stanford.edu/cgi-bin/wa?SUBED1=LCD-CVS&A=1