Author: mgraham Date: Thu Oct 30 08:16:13 2014 New Revision: 1342 Log: Straight track fitting Added: java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTrack2DFitter.java java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTrackConfirmerExtender.java java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTrackFitter.java java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTrackReconDriver.java - copied, changed from r1321, java/trunk/tracking/src/main/java/org/hps/recon/tracking/TrackerReconDriver.java java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTracker.java Removed: java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/HitOnTrackChecker.java java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTrack.java java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/TrackChecker.java java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/TrackCollectionUtilities.java Modified: java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTrackFinder.java Added: java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTrack2DFitter.java ============================================================================= --- java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTrack2DFitter.java (added) +++ java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTrack2DFitter.java Thu Oct 30 08:16:13 2014 @@ -0,0 +1,131 @@ +package org.hps.recon.tracking.nobfield; + +import hep.physics.matrix.SymmetricMatrix; +import java.util.HashMap; +import java.util.List; +import java.util.Map; +import org.lcsim.fit.helicaltrack.HelicalTrackFit; +import org.lcsim.fit.helicaltrack.HelicalTrackHit; +import org.lcsim.fit.helicaltrack.HitUtils; +import org.lcsim.fit.helicaltrack.MultipleScatter; +import org.lcsim.fit.line.SlopeInterceptLineFit; +import org.lcsim.fit.line.SlopeInterceptLineFitter; + +/** + * + * @author mgraham + */ +public class StraightTrack2DFitter { + + SlopeInterceptLineFitter _lfitter = new SlopeInterceptLineFitter(); + HelicalTrackFit _fit; + + /** + * Status of the HelicalTrackFit. + */ + public enum FitStatus { + + /** + * Successful Fit. + */ + Success, + /** + * Inconsistent seed hits + */ + InconsistentSeed, + /** + * s-z line fit failed. + */ + SZLineFitFailed, + /** + * ZSegmentFit failed. + */ + XYLineFitFailed + }; + + public void StraightTrack2DFitter() { + + } + + public FitStatus fit(List<HelicalTrackHit> hits) { + Map<HelicalTrackHit, MultipleScatter> msmap = new HashMap<HelicalTrackHit, MultipleScatter>(); + return fit(hits, msmap, null); + } + + public FitStatus fit(List<HelicalTrackHit> hits, Map<HelicalTrackHit, MultipleScatter> msmap, HelicalTrackFit oldfit) { + // Check that we have at least 3 hits + boolean success = false; + int nhit = hits.size(); + if (nhit < 3) + return FitStatus.InconsistentSeed; + + // Create the objects that will hold the fit output + double[] chisq = new double[2]; + int[] ndof = new int[2]; + double[] par = new double[5]; + SymmetricMatrix cov = new SymmetricMatrix(5); + + // Setup for the line fit + double[] s = new double[nhit]; + double[] z = new double[nhit]; + double[] y = new double[nhit]; + double[] x = new double[nhit]; + double[] dz = new double[nhit]; + double[] dy = new double[nhit]; + Map<HelicalTrackHit, Double> smap = new HashMap<>(); + + // Store the coordinates and errors for the XY line fit + for (int i = 0; i < nhit; i++) { + HelicalTrackHit hit = hits.get(i); + y[i] = hit.y(); + dy[i] = HitUtils.zres(hit, msmap, oldfit); + double drphi_ms = 0; + if (msmap.containsKey(hit)) + drphi_ms = msmap.get(hit).drphi(); + double dyHitSq = hit.getCorrectedCovMatrix().e(1, 1); + dy[i] = Math.sqrt(dyHitSq + drphi_ms * drphi_ms); + x[i] = hit.x(); + } + // Call the line fitter and check for success + success = _lfitter.fit(x, y, dy, nhit); + if (!success) + return FitStatus.XYLineFitFailed; + SlopeInterceptLineFit xyFit = _lfitter.getFit(); + par[0] = xyFit.intercept(); + par[1] = xyFit.slope(); + cov.setElement(0, 0, Math.pow(xyFit.interceptUncertainty(), 2)); + cov.setElement(1, 1, Math.pow(xyFit.slopeUncertainty(), 2)); + cov.setElement(0, 1, Math.pow(xyFit.covariance(), 2)); + chisq[0] = xyFit.chisquared(); + ndof[0] = xyFit.ndf(); + // Store the coordinates and errors for the SZ line fit + for (int i = 0; i < nhit; i++) { + HelicalTrackHit hit = hits.get(i); + z[i] = hit.z(); + dz[i] = HitUtils.zres(hit, msmap, oldfit); //MG this works even in msmap is null + s[i] = Math.sqrt(hit.x() * hit.x() + hit.y() * hit.y()); + smap.put(hit, s[i]); + } + // Call the line fitter and check for success + success = _lfitter.fit(s, z, dz, nhit); + if (!success) + return FitStatus.SZLineFitFailed; + SlopeInterceptLineFit szFit = _lfitter.getFit(); + par[3] = szFit.intercept(); + par[4] = szFit.slope(); + cov.setElement(3, 3, Math.pow(szFit.interceptUncertainty(), 2)); + cov.setElement(4, 4, Math.pow(szFit.slopeUncertainty(), 2)); + cov.setElement(3, 4, Math.pow(szFit.covariance(), 2)); + chisq[1] = szFit.chisquared(); + ndof[1] = szFit.ndf(); + par[2] = 1e-50; + // Create the HelicalTrackFit for this helix + _fit = new HelicalTrackFit(par, cov, chisq, ndof, smap, msmap); + return FitStatus.Success; + } + + public HelicalTrackFit getFit() { + return _fit; + } + +} Added: java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTrackConfirmerExtender.java ============================================================================= --- java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTrackConfirmerExtender.java (added) +++ java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTrackConfirmerExtender.java Thu Oct 30 08:16:13 2014 @@ -0,0 +1,356 @@ +/* + * StraightTrackConfirmerExtender.java + */ +package org.hps.recon.tracking.nobfield; + +import java.util.ArrayList; +import java.util.Collections; +import java.util.HashMap; +import java.util.LinkedList; +import java.util.List; +import java.util.Map; +import org.lcsim.fit.helicaltrack.HelicalTrackFit; +import org.lcsim.fit.helicaltrack.HelicalTrackHit; +import org.lcsim.recon.tracking.seedtracker.HitManager; +import org.lcsim.recon.tracking.seedtracker.MergeSeedLists; +import org.lcsim.recon.tracking.seedtracker.SeedCandidate; +import org.lcsim.recon.tracking.seedtracker.SeedLayer; +import org.lcsim.recon.tracking.seedtracker.SeedStrategy; +import org.lcsim.recon.tracking.seedtracker.SortHits; +import org.lcsim.recon.tracking.seedtracker.SortLayers; +import org.lcsim.recon.tracking.seedtracker.diagnostic.ISeedTrackerDiagnostics; + +/** + * The StraightTrackConfirmerExtender class attempts to add hits to an input seed. While the + * algorithms used for the confirm and extend phases are identical, there are + * small differences in procedures when there are no more tracker layers to check. + * The confirm phase simply outputs a list of SeedCandidates that have at least + * the minimum number of confirm hits added to the track. The extend phases + * completes track finding and imposes the merge criteria to eliminate inferior + * track candidates when a pair of candidates shares more than one hit. + * + * @author cozzy, Richard Partridge + * Modified for HPS straight track fitting by Matt Graham 10/29/2014 + */ +public class StraightTrackConfirmerExtender { + + public enum Task { + + CONFIRM, + EXTEND; + } + private int _nfit; + private int _maxfit = 1000000000; // initialize maximum number of fits to 10^9 + private HitManager _hmanager; + private StraightTrackFitter _fitter; + private MergeSeedLists _merger; + private List<SeedCandidate> _result; + private ISeedTrackerDiagnostics _diag = null; + + + /** + * Constructor for the ConfirmerExtender class. + * + * @param hitmanager hit manager that provides access to hits sorted by layer / sector + * @param helixfitter helix fitter + */ + public StraightTrackConfirmerExtender(HitManager hitmanager, StraightTrackFitter helixfitter) { + + _hmanager = hitmanager; + _fitter = helixfitter; + _merger = new MergeSeedLists(); + } + + /** + * Retrieves the result of the confirmation or extension. + * + * @return result The list of confirmed/extended seeds.. + */ + public List<SeedCandidate> getResult() { + + return _result; + } + /** + * Set the maximum number of fit trials for a given seed to be confirmed/extended. + * + * @param maxfit maximum number of trials + */ + public void setMaxFit(int maxfit) { + _maxfit = maxfit; + } + + /** + * Get the number of fit trials for the last confirm/extend. + * + * @return number of helix fits + */ + public int getNFit() { + return _nfit; + } + + /** + * Set the diagnostic class to be used (null for no diagnostics). + * + * @param d diagnostic class + */ + public void setDiagnostics(ISeedTrackerDiagnostics d) { + _diag = d; + _merger.setDiagnostic(_diag); + } + + /** + * Try to confirm a seed using a specified strategy. The strategy specifies + * the layers to use in trying to confirm the seed as well as the minimum number + * of confirm layers that were added to the seed. Typically, there will be a + * single confirm layer that is required for a seed to be confirmed. Confirmed + * seeds must also pass the chisq cut and other constraints specified in the strategy. + * + * @param seed seed to be confirmed + * @param strategy strategy to use for confirming this seed + * @param bfield magnetic field + * @return true if seed was confirmed + */ + public boolean Confirm(SeedCandidate seed, SeedStrategy strategy, double bfield) { + + // Create a list to hold the confirmed / extended seeds + _result = new ArrayList<SeedCandidate>(); + + // Establish which layers are to be checked + seed.setUncheckedLayers(strategy.getLayers(SeedLayer.SeedType.Confirm)); + + // Process the seed + doTask(seed, Task.CONFIRM, strategy, bfield); + + // Return true if we found at least one confirming seed candidate + return _result.size() > 0; + } + + /** + * Try to extend a seed using a specified strategy. The strategy specifies + * the layers to use in extending the seed as well as the minimum number of + * layers required for a seed to be a track candidate. Any track candidates + * found as a result of the extend operation are merged with the list of + * track candidates found so far to eliminate inferior fits when a pair of + * track candidates shares more than one hit. + * + * @param seed seed to be extended + * @param strategy strategy to use + * @param bfield magnetic field + * @param foundseeds list of track candidates found so far + */ + public void Extend(SeedCandidate seed, SeedStrategy strategy, double bfield, + List<SeedCandidate> foundseeds) { + + // Initialize the list of extended seed candidates to those already found + _result = foundseeds; + + // Establish which layers are to be checked + seed.setUncheckedLayers(strategy.getLayers(SeedLayer.SeedType.Extend)); + + // Extend the seed and return + doTask(seed, Task.EXTEND, strategy, bfield); + return; + } + + /** + * Perform the confirm or extend step. + * + * @param inputseed seed to be confirmed/extended + * @param task confirm or extend (enum) + * @param strategy strategy to use + * @param bfield magnetic field + */ + private void doTask(SeedCandidate inputseed, Task task, SeedStrategy strategy, double bfield) { + + // Initialize the counter for the number of fits performed on this seed + _nfit = 0; + +// // Instantiate the fast hit checker +// FastCheck checker = new FastCheck(strategy, bfield, _diag); +// if(this._applySectorBinning) checker.setDoSectorBinCheck(this._hmanager.getSectorManager()); + + // Calculate the minimum number of hits to succeed, retrieve the chisq cuts + int minhits = strategy.getMinHits(); + if (task == Task.CONFIRM) minhits = strategy.getMinConfirm() + 3; + double badhitchisq = strategy.getBadHitChisq(); + double maxchisq = strategy.getMaxChisq(); + + // Create a LIFO queue of seeds to be searched for a confirmation/extension + // hit (note that a LIFO queue is used to minimize memory usage) + LinkedList<SeedCandidate> seedlist = new LinkedList<SeedCandidate>(); + + // The bestseed is a SeedCandidate the meets the requirements for becoming + // a track, shares at least one hit with the inputseed, and has been deemed + // the best such seed by the track merging criteria + // + // Initialize the best seed to null + SeedCandidate bestseed = null; + + // If we have already found track candidates, check for duplicates + // that share hits with the seed, finding the best such duplicate candidate. + for (SeedCandidate trkcand : _result) { + if (_merger.isDuplicate(inputseed, trkcand)) + bestseed = findBestCandidate(trkcand, bestseed, strategy); + } + + // Create a map between the SeedLayers to be checked and a list of hits on the layer to check + Map<SeedLayer, List<HelicalTrackHit>> hitmap = new HashMap<SeedLayer, List<HelicalTrackHit>>(); + + // Loop over the layers to be checked + for (SeedLayer lyr : inputseed.getUncheckedLayers()) { + + hitmap.put(lyr, _hmanager.getTrackerHits(lyr)); + } + + // Create a list of layers that have hits to check + List<SeedLayer> lyrlist = new ArrayList<SeedLayer>(); + lyrlist.addAll(hitmap.keySet()); + + // Sort the layers in order of increasing number of hits + SortLayers lyrsort = new SortLayers(hitmap); + Collections.sort(lyrlist, lyrsort); + + // Store the layers to check in the seed + inputseed.setUncheckedLayers(lyrlist); + + // Start with the input seed + seedlist.add(inputseed); + + // Keep looping until we have fully processed all seed candidates + while (seedlist.size() > 0) { + + // If we have exceeded the maximum number of fits, print warning and stop processing seed candidates + if (_nfit > _maxfit) { + System.out.println("Maximum number of fits exceeded in "+task.toString()+" step"); + if (bestseed == null) { + System.out.println("No track candidates are associated with the seed hits"); + } else { + System.out.println("Track candidate with "+bestseed.getHits().size()+ + " hits and chisq of "+bestseed.getHelix().chisqtot()+ + " associated with the seed hits"); + } + break; + } + + // Pull the last seed off the queue (use a LIFO queue to minimize queue length) + SeedCandidate seed = seedlist.poll(); + + // Check if there are enough unchecked layers to meet the minimum number of hits + int lyrsleft = seed.getUncheckedLayers().size(); + int possiblehits = lyrsleft + seed.getHits().size(); + if (possiblehits < minhits) continue; + + // If there is a best fit candidate, see if there is still a chance of beating it + if (bestseed != null) { + + // If the maximimum hits we can achieve is >1 fewer than the best fit, skip this candidate + int besthits = bestseed.getHits().size(); + if (possiblehits < besthits - 1) continue; + + // If the maximum hits we can achieve equals the best fit, skip if we have a worse chi2 + double chisq = seed.getHelix().chisqtot(); + double bestchisq = seed.getHelix().chisqtot(); + if ((possiblehits == besthits) && chisq > bestchisq) continue; + + // If the maximum hits we can achieve is 1 fewer than the best fit, skip if the bad hit criteria can't be met + if ((possiblehits == besthits - 1) && (chisq > bestchisq - badhitchisq)) continue; + } + + // See if there are any layers left for confirm/extend + if (lyrsleft == 0) { + + // Take final action on this seed + if (task == Task.CONFIRM) { + + // No more layers and min hit requirement is met, seed is confirmed + _result.add(seed); + + } else if (task == Task.EXTEND) { + + // Merge the seed into the list of extended seeds + boolean merged = _merger.Merge(_result, seed, strategy); + + // If the seed survived the merge, make it our new best candidate + if (merged) bestseed = findBestCandidate(seed, bestseed, strategy); + } + + // Done with this seed + continue; + } + + // Pull the next layer off the queue + SeedLayer lyr = seed.getNextLayer(); + HelicalTrackFit helix = seed.getHelix(); + + // Retrieve the chisq for the last fit and initialize the best fit chisq for this layer + double oldchisq = helix.chisqtot(); + double oldcirclechisq = helix.chisq()[0]; + double chisqbest = 1.e99; + + // Get the list of hits to check for this layer and sort them by x-y distance from current helix + List<HelicalTrackHit> hitlist = hitmap.get(lyr); + SortHits comp = new SortHits(helix); + Collections.sort(hitlist, comp); + + // Loop over the sorted hits in this layer + for (HelicalTrackHit hit : hitlist) { + + // Make a test seed including the new hit + SeedCandidate test = new SeedCandidate(seed); + test.addHit(hit); + +// // Check that this hit is potentially viable +// if (!checker.CheckHitSeed(hit, seed)) { +// if (_diag != null) _diag.fireCheckHitFailed(hit, test); +// continue; +// } + + // Fit the test seed + boolean success = _fitter.FitCandidate(test, strategy); + _nfit++; + + + // Check if the fit was successful + if (success) { + + // Success - attach the fit to the test seed + HelicalTrackFit newhelix = _fitter.getHelix(); + test.setHelix(newhelix); + + // Add the seed to the LIFO queue of seed candidates and update the best chisq + seedlist.addFirst(test); + chisqbest = Math.min(chisqbest, newhelix.chisqtot()); + + } + } + + // Finished checking hits in the current layer. If all the fit trials for + // this layer are potentially bad hits, include the starting seed (less the + // current layer, which was popped off the layer queue) in the seed list. + if (chisqbest - oldchisq > strategy.getBadHitChisq()) seedlist.addFirst(seed); + } + + // Finished looping over the seeds in the LIFO candidate queue - we are done! + return; + } + + /** + * Check two track candidates and return the best one subject to the merging criteria. + * + * @param trial new candidate to try + * @param oldbest previous best candidate (or null if no best candidate has been found) + * @param strategy strategy in use + * @return best track candidate + */ + private SeedCandidate findBestCandidate(SeedCandidate trial, SeedCandidate oldbest, SeedStrategy strategy) { + + // If the old best candidate is null, return the trial candidate + if (oldbest == null) return trial; + + // If the trial candidate is better, return it + if (_merger.isBetter(trial, oldbest, strategy)) return trial; + + // If no improvement, return the old candidate + return oldbest; + } +} Modified: java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTrackFinder.java ============================================================================= --- java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTrackFinder.java (original) +++ java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTrackFinder.java Thu Oct 30 08:16:13 2014 @@ -1,272 +1,252 @@ +/* + * SeedTrackFinder.java + * + * Created on January 22, 2008, 9:39 AM + * + */ package org.hps.recon.tracking.nobfield; import java.util.ArrayList; +import java.util.HashSet; import java.util.List; -import org.hps.recon.tracking.HitCollectionUtilites; -import org.hps.recon.tracking.nobfield.TrackCollectionUtilities; -import org.lcsim.event.EventHeader; -import org.lcsim.event.Track; -import org.lcsim.event.TrackerHit; +import java.util.Set; +import org.lcsim.event.MCParticle; import org.lcsim.fit.helicaltrack.HelicalTrackHit; -import org.lcsim.fit.line.SlopeInterceptLineFit; -import org.lcsim.fit.line.SlopeInterceptLineFitter; -import org.lcsim.geometry.Detector; -import org.lcsim.util.Driver; +import org.lcsim.recon.tracking.seedtracker.HitManager; +import org.lcsim.recon.tracking.seedtracker.SeedCandidate; +import org.lcsim.recon.tracking.seedtracker.SeedLayer; +import org.lcsim.recon.tracking.seedtracker.SeedStrategy; +import org.lcsim.recon.tracking.seedtracker.TrackCheck; +import org.lcsim.recon.tracking.seedtracker.diagnostic.ISeedTrackerDiagnostics; /** * - * @author mgraham + * @author Richard Partridge + * @version 1.0 */ -public class StraightTrackFinder extends Driver { - - // Debug flag. - private boolean debug = true; - // Tracks found across all events. - int ntracks = 0; - // Number of events processed. - int nevents = 0; - // Cache detector object. - Detector detector = null; - // Tracking strategies resource path. - private String strategyResource = "HPS-Test-4pt1.xml"; - // Output track collection. - private String trackCollectionName = "StraightTracks"; - // HelicalTrackHit input collection. - private String stInputCollectionName = "HelicalTrackHits"; - // Include MS (done by removing XPlanes from the material manager results) - private boolean includeMS = true; - // number of repetitive fits on confirmed/extended tracks - private int _iterativeConfirmed = 3; - // use HPS implementation of material manager - private boolean _useHPSMaterialManager = true; - - private TrackChecker checkerTrack = new TrackChecker(); - private HitOnTrackChecker checkerHOT = new HitOnTrackChecker(); - - private SlopeInterceptLineFitter _lfitter = new SlopeInterceptLineFitter(); - - public void setDebug(boolean debug) { - this.debug = debug; - } - - /** - * Set the tracking strategy resource. - * - * @param strategyResource The absolute path to the strategy resource in the - * hps-java jar. - */ - public void setStrategyResource(String strategyResource) { - this.strategyResource = strategyResource; - } - - public void setInputHitCollectionName(String inputHitCollectionName) { - this.stInputCollectionName = inputHitCollectionName; - } - - public void setTrackCollectionName(String trackCollectionName) { - this.trackCollectionName = trackCollectionName; - } - - public void setIncludeMS(boolean incMS) { - this.includeMS = incMS; - } - - /** - * Set to enable the use of the HPS material manager implementation - * - * @param useHPSMaterialManager switch - */ - public void setUseHPSMaterialManager(boolean useHPSMaterialManager) { - this._useHPSMaterialManager = useHPSMaterialManager; - } - - @Override - public void detectorChanged(Detector detector) { - // Cache Detector object. - this.detector = detector; -// initialize(); - super.detectorChanged(detector); - } - - @Override - public void process(EventHeader event) { - if (!event.hasCollection(HelicalTrackHit.class, stInputCollectionName)) - return; - - List<HelicalTrackHit> allHits = event.get(HelicalTrackHit.class, stInputCollectionName); - - List<List<HelicalTrackHit>> splitTopBot = HitCollectionUtilites.SplitTopBottomHits(allHits); - // will always have top(=0) and bottom(=1) lists (though they may be empty) - List<HelicalTrackHit> topHits = splitTopBot.get(0); - List<HelicalTrackHit> bottomHits = splitTopBot.get(1); - //a simple strategy...eventually implement SeedTracker strategies? - int nTotLayers = 6; - int nSeed = 3; - int nExtra = nTotLayers - nSeed; - int[] seedStrategy = {1, 3, 5}; - int[] extendStrategy = {7, 9, 11}; - int minHits = 4; - - TrackChecker checkerTrack = new TrackChecker(); - HitOnTrackChecker checkerHOT = new HitOnTrackChecker(); - -// List<StraightTrack> seeds = getSeeds(seedStrategy, topHits); - List<StraightTrack> seeds = getSeeds(seedStrategy, allHits); - System.out.println("Found " + seeds.size() + " seeds"); - - List<StraightTrack> extendedSeeds = new ArrayList<>(); - for (StraightTrack seed : seeds) - extendTrack(extendStrategy, 0, seed, allHits, extendedSeeds); -// extendTrack(extendStrategy, 0, seed, topHits, extendedSeeds); - - System.out.println("Prepruning :Found " + extendedSeeds.size() + " extended seeds"); - - //remove tracks with more than m overlaping hits...pick best chi2 - //... - List<StraightTrack> finalTracks = new ArrayList<>(); - for (StraightTrack track : extendedSeeds) { - boolean isbest = TrackCollectionUtilities.pruneTrackList((ArrayList<Track>) (ArrayList) extendedSeeds, track, 1); - if (isbest) - finalTracks.add(track); - } - - System.out.println("Postpruning :Found " + finalTracks.size() + " extended seeds"); - event.put(trackCollectionName, finalTracks); - } - - public SlopeInterceptLineFit FitToLine(List<HelicalTrackHit> hits, int projection) { - SlopeInterceptLineFit _lfit; - int npix = hits.size(); - double[] s = new double[npix]; - double[] q = new double[npix]; - double[] dq = new double[npix]; - - // Store the coordinates and errors for the line fit - for (int i = 0; i < npix; i++) { - HelicalTrackHit hit = hits.get(i); - s[i] = hit.z();//probably isn't quite right...track length is not z - if (projection == 0) { //do x vs z; - q[i] = hit.x(); - dq[i] = Math.sqrt(hit.getCorrectedCovMatrix().e(0, 0)); - } else { - q[i] = hit.y(); - dq[i] = Math.sqrt(hit.getCorrectedCovMatrix().e(1, 1)); - } - } - - // Call the line fitter and check for success - boolean success = _lfitter.fit(s, q, dq, npix); - - if (!success) - System.out.println("Something is broken in the line fit"); - // Save the line fit, chi^2, and DOF - _lfit = _lfitter.getFit(); - return _lfit; - - } - - private StraightTrack makeTrack(List<HelicalTrackHit> hits, SlopeInterceptLineFit xfit, SlopeInterceptLineFit yfit) { - StraightTrack track = new StraightTrack(); - double[] pars = {-99, -99, -99, -99, -99};//this needs to have 5 fields to implement Track - pars[0] = xfit.intercept(); - pars[1] = xfit.slope(); - pars[2] = yfit.intercept(); - pars[3] = yfit.slope(); - track.setTrackParameters(pars); - track.setChi2(xfit.chisquared(), yfit.chisquared()); - track.setNDF(xfit.ndf()+yfit.ndf()); - for (TrackerHit hit : hits) - track.addHit(hit); - // TODO: set convariance, - return track; - } - - private StraightTrack makeTrack(List<HelicalTrackHit> hits) { - SlopeInterceptLineFit xfit = FitToLine(hits, 0); - SlopeInterceptLineFit yfit = FitToLine(hits, 1); - if (debug) - System.out.println("xfit = " + xfit.toString()); - if (debug) - System.out.println("yfit = " + yfit.toString()); - return makeTrack(hits, xfit, yfit); - } - - /* - * Get all seed combinations that make sense (pass checkSeed) - * currently, just assume there are 3 seed layers (don't have to be first 3 though. - */ - private List<StraightTrack> getSeeds(int[] seedLayers, List<HelicalTrackHit> hits) { - List<StraightTrack> seeds = new ArrayList<>(); - int nseeds = seedLayers.length; - if (nseeds == 3)//TODO ... set this up so that this works for arbitrary nseeds...use recursion - for (HelicalTrackHit h1 : HitCollectionUtilites.GetSortedHits(hits, seedLayers[0])) { - if (debug) - System.out.println(h1.toString()); - for (HelicalTrackHit h2 : HitCollectionUtilites.GetSortedHits(hits, seedLayers[1])) { - if (debug) - System.out.println(h2.toString()); - for (HelicalTrackHit h3 : HitCollectionUtilites.GetSortedHits(hits, seedLayers[2])) { - if (debug) - System.out.println(h3.toString()); - //make a 3-hit test track...see if it passes CheckTrack - List<HelicalTrackHit> testTrack = new ArrayList<HelicalTrackHit>(); - testTrack.add(h1); - testTrack.add(h2); - testTrack.add(h3); - StraightTrack trk = makeTrack(testTrack); - if (!checkerTrack.checkSeed(trk)) - break; - seeds.add(trk); +public class StraightTrackFinder { + + private HitManager _hitmanager; + private StraightTrackFitter _helixfitter; + private StraightTrackConfirmerExtender _confirmer; + private List<SeedCandidate> _trackseeds; + private ISeedTrackerDiagnostics _diag = null; + private Set<MCParticle> _seededmcp; + private Set<MCParticle> _confirmedmcp; + TrackCheck _trackCheck; // set by SeedTracker + private boolean _debug = false; + + + /** + * Creates a new instance of SeedTrackFinder + */ + public StraightTrackFinder(HitManager hitmanager, StraightTrackFitter helixfitter) { + + // Save the pointers to the hit manager and helix fitter classes + _hitmanager = hitmanager; + _helixfitter = helixfitter; + + // Instantiate the Confirmer/Extender and Seed Candidate merging classes + _confirmer = new StraightTrackConfirmerExtender(_hitmanager, _helixfitter); + + // Create a list of track seeds that have been found + _trackseeds = new ArrayList<SeedCandidate>(); + + // Create a set of MC Particles that have been seeded, confirmed + _seededmcp = new HashSet<MCParticle>(); + _confirmedmcp = new HashSet<MCParticle>(); + + } + + public void setDiagnostic(ISeedTrackerDiagnostics d) { + + // Setup the diagnostics for this class and the classes used by this class + _diag = d; + _confirmer.setDiagnostics(_diag); + } + + public boolean FindTracks(SeedStrategy strategy, double bfield) { + + // Instantiate the fast hit checker + // FastCheck checker = new FastCheck(strategy, bfield, _diag); + // if(_applySectorBinning) checker.setDoSectorBinCheck(_hitmanager.getSectorManager()); +// SeedSectoring ss = new SeedSectoring(_hitmanager, strategy, bfield, _applySectorBinning); + List<SeedLayer> seeds = strategy.getLayers(SeedLayer.SeedType.Seed); + +// List<List<Sector>> sslist = ss.SeedSectors(); +// if(_debug) +// System.out.println(this.getClass().getSimpleName()+": number of SeedSectors="+sslist.size()); + // Loop over the valid sector combinations + // for (List<Sector> slist : sslist) { + // Loop over the first seed layer + for (HelicalTrackHit hit1 : _hitmanager.getTrackerHits(seeds.get(0))) + + // Loop over the second seed layer and check that we have a hit pair consistent with our strategy + for (HelicalTrackHit hit2 : _hitmanager.getTrackerHits(seeds.get(1))) { + + // Call _trackCheck if set + if (_trackCheck != null) { + SeedCandidate tempseed = new SeedCandidate(strategy, bfield); + tempseed.addHit(hit1); + tempseed.addHit(hit2); + if (!_trackCheck.checkSeed(tempseed)) + continue; + } + +// // Check if the pair of hits is consistent with the current strategy +// if (!checker.TwoPointCircleCheck(hit1, hit2, null)) { +// if (_diag != null) _diag.fireCheckHitPairFailed(hit1, hit2); +// continue; +// } + // Loop over the third seed layer and check that we have a hit triplet consistent with our strategy + for (HelicalTrackHit hit3 : _hitmanager.getTrackerHits(seeds.get(2))) { + // Call _trackCheck if set + if (_trackCheck != null) { + SeedCandidate tempseed2 = new SeedCandidate(strategy, bfield); + tempseed2.addHit(hit1); + tempseed2.addHit(hit3); + if (!_trackCheck.checkSeed(tempseed2)) + continue; + + SeedCandidate tempseed3 = new SeedCandidate(strategy, bfield); + tempseed3.addHit(hit2); + tempseed3.addHit(hit3); + if (!_trackCheck.checkSeed(tempseed3)) + continue; } + + // Form a seed candidate from the seed hits + SeedCandidate seed = new SeedCandidate(strategy, bfield); + seed.addHit(hit1); + seed.addHit(hit2); + seed.addHit(hit3); + +// // Check if the triplet of hits is consistent with the current strategy +// if (!checker.ThreePointHelixCheck(hit1, hit2, hit3)) { +// +// if (_diag != null) { +// if (seed.isTrueSeed()) +// _diag.fireCheckHitTripletFailed(hit1, hit2, hit3); +// } +// continue; +// } + // Form a seed candidate from the seed hits + // If it's a true seed, add the MC Particle to those that were seeded + if (_diag != null) + if (seed.isTrueSeed()) + _seededmcp.addAll(seed.getMCParticles()); + + if (_debug) + System.out.println(this.getClass().getSimpleName() + ": fit the candidate"); + + // See if we can fit a helix to this seed candidate + boolean success = _helixfitter.FitCandidate(seed, strategy); + + if (!success) + continue; + + if (_debug) + System.out.println(this.getClass().getSimpleName() + ": fit success"); + + // Save the helix fit + seed.setHelix(_helixfitter.getHelix()); + + // Check the seed - hook for plugging in external constraint + if (_trackCheck != null) + if (!_trackCheck.checkSeed(seed)) + continue; + + // See if we can confirm this seed candidate + success = _confirmer.Confirm(seed, strategy, bfield); + if (!success) + continue; + + if (_debug) + System.out.println(this.getClass().getSimpleName() + ": confirmed seed"); + + // Confirmed a seed - if it's a true seed, add the MC Particle to those that were confirmed + if (_diag != null) + if (seed.isTrueSeed()) + _confirmedmcp.addAll(seed.getMCParticles()); + + if (_debug) + System.out.println(this.getClass().getSimpleName() + ": try to extend"); + + // Try to extend each confirmed seed candidates to make a track candidate + List<SeedCandidate> confirmedlist = _confirmer.getResult(); + for (SeedCandidate confirmedseed : confirmedlist) + + // See if we can extend this seed candidate + _confirmer.Extend(confirmedseed, strategy, bfield, _trackseeds); } } - return seeds; - } - - /* - * recursively extend the seeds through all of the extend layers.. - * ...I think this should work... - */ - private void extendTrack(int[] extendLayers, int n, StraightTrack origTrack, List<HelicalTrackHit> hits, List<StraightTrack> trackList) { - if (n >= extendLayers.length) { - if (debug) - System.out.println("Done finding this track through all " + n + " extra layers"); - trackList.add(origTrack); - return; - } - - boolean cannotExtendThisLayer = true; - if (debug) - System.out.println("Extending to layer " + extendLayers[n]); - for (HelicalTrackHit h : HitCollectionUtilites.GetSortedHits(hits, extendLayers[n])) { - //let's see if this hit makes sense to add to original track - if (!checkerHOT.checkNewHit(origTrack, h)) - continue; - - List<TrackerHit> origHits = origTrack.getTrackerHits(); - //make a new list and cast them as HelicalTrackHits (Track only stores TrackerHits) - List<HelicalTrackHit> newHits = new ArrayList<>(); - for (TrackerHit oh : origHits) { - HelicalTrackHit hoh = (HelicalTrackHit) oh; - System.out.println(hoh.getPosition()[0]); - newHits.add(hoh); - } - //add the new hit to the list & make new track - newHits.add(h); - StraightTrack newTrack = makeTrack(newHits); - //check the new track after we've added this hit - if (!checkerTrack.checkTrack(newTrack)) - continue; - cannotExtendThisLayer = false; - //extend again to the next layer - extendTrack(extendLayers, n + 1, newTrack, hits, trackList); - } - - //didn't find any hits in this layer that match the track...but let's try the next one - if (cannotExtendThisLayer) - extendTrack(extendLayers, n + 1, origTrack, hits, trackList); - - return; + +// // Done with track finding for this strategy +// if (_diag != null) +// _diag.fireFinderDone(_trackseeds, _seededmcp); + return _trackseeds.size() > 0; + } + + /** + * Return the list of track candidates. + * + * @return track candidates + */ + public List<SeedCandidate> getTrackSeeds() { + return _trackseeds; + } + + /** + * Clear the list of track candidates accumulated from previous calls to + * SeedTrackFinder (typically done before starting a new event). + */ + public void clearTrackSeedList() { + _trackseeds.clear(); + _seededmcp.clear(); + _confirmedmcp.clear(); + } + + /** + * Set the maximum number of fits for a given seed in a confirm or extend + * step. + * + * @param maxfits maximum number of fits + */ + public void setMaxFit(int maxfits) { + _confirmer.setMaxFit(maxfits); + } + + /** + * Return the list of MCParticles that formed valid 3-hit seeds. + * + * @return list of seeded MCParticles + */ + public Set<MCParticle> getSeededMCParticles() { + return _seededmcp; + } + + /** + * Return the list of confirmed MCParticles. + * + * @return confirmed MCParticles + */ + public Set<MCParticle> getConfirmedMCParticles() { + return _confirmedmcp; + } + + /** + * Return the StraightTrackConfirmerExtender + * + * @return confirmer/extender object + * + */ + public StraightTrackConfirmerExtender getConfirmer() { + return _confirmer; + } + + public void setDebug(boolean debug) { + System.out.println("Setting " + this.getClass().getSimpleName() + " debug to " + debug); + _debug = debug; } } Added: java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTrackFitter.java ============================================================================= --- java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTrackFitter.java (added) +++ java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTrackFitter.java Thu Oct 30 08:16:13 2014 @@ -0,0 +1,261 @@ +/* + * HelixFitter.java + * + * Created on January 22, 2008, 9:25 AM + * + */ +package org.hps.recon.tracking.nobfield; + +import java.util.List; +import java.util.Map; +import org.hps.recon.tracking.MaterialManager; +import org.hps.recon.tracking.MultipleScattering; +import org.hps.recon.tracking.nobfield.StraightTrack2DFitter.FitStatus; +import org.lcsim.fit.circle.CircleFit; +import org.lcsim.fit.helicaltrack.HelicalTrackCross; +import org.lcsim.fit.helicaltrack.HelicalTrackFit; +import org.lcsim.fit.helicaltrack.HelicalTrackHit; +import org.lcsim.fit.helicaltrack.MultipleScatter; +import org.lcsim.fit.line.SlopeInterceptLineFit; +import org.lcsim.fit.line.SlopeInterceptLineFitter; +import org.lcsim.fit.zsegment.ZSegmentFit; +import org.lcsim.recon.tracking.seedtracker.ConstrainHelix; +import org.lcsim.recon.tracking.seedtracker.SeedCandidate; +import org.lcsim.recon.tracking.seedtracker.SeedStrategy; +import org.lcsim.recon.tracking.seedtracker.TrackCheck; +import org.lcsim.recon.tracking.seedtracker.diagnostic.ISeedTrackerDiagnostics; + +/** + * + * @author Richard Partridge + * @version 1.0 + */ +public class StraightTrackFitter { + + private StraightTrack2DFitter _fitter = new StraightTrack2DFitter(); + protected MultipleScattering _scattering; + private HelicalTrackFit _helix;//still use HelicalTrackFit... + private MaterialManager _materialmanager; + //private ConstrainHelix _constrain; + private double _bfield = 0.; + private CircleFit _circlefit; + private SlopeInterceptLineFit _linefit; + private ZSegmentFit _zsegmentfit; + private FitStatus _status; + private ISeedTrackerDiagnostics _diag = null; + TrackCheck _trackCheck; // set by SeedTracker + private boolean _debug = false; + + /** + * Creates a new instance of StraightTrackFitter + */ + public StraightTrackFitter(MaterialManager materialmanager) { + _materialmanager = materialmanager; + + _scattering = new MultipleScattering(_materialmanager); + + // _constrain = new ConstrainHelix(); + } + + // public boolean FitCandidate(SeedCandidate seed, SeedStrategy strategy) { + public boolean FitCandidate(SeedCandidate seed, SeedStrategy strategy) {// change this to something reasonable + if (_debug) + System.out.printf("%s: FitCandidate() with strategy: \"%s\"\n", this.getClass().getSimpleName(), strategy.getName()); + + // Initialize fit results to null objects + _helix = null; + + // Check that we have set the magnetic field + if (_bfield != 0.) + throw new RuntimeException("B Field should be zero for the straight track fitter"); + + // Set the tolerance in the fitter + // _fitter.setTolerance(Math.sqrt(strategy.getMaxChisq())); + // Retrieve list of hits to be fit + List<HelicalTrackHit> hitlist = seed.getHits(); + + if (_debug) { + System.out.println(this.getClass().getSimpleName() + ": hitlist size " + hitlist.size() + ":"); + double z_prev = -99999999.9; + for (HelicalTrackHit hit : hitlist) { + System.out.printf("%s: (%.2f,%.2f,%.2f) corrected %s\n", + this.getClass().getSimpleName(), hit.getPosition()[0], hit.getPosition()[1], hit.getPosition()[2], + hit.getCorrectedPosition().toString()); + if (z_prev < -99999999.0) + z_prev = hit.getPosition()[2]; + else { + if (Math.signum(z_prev) != Math.signum(hit.getPosition()[2])) + System.out.printf("%s: topBotHits in event\n", this.getClass().getSimpleName()); + z_prev = hit.getPosition()[2]; + } + } + } + // Retrieve the old helix + HelicalTrackFit oldhelix = seed.getHelix(); + + // If this is the candidate's first helix fit, first do a fit without MS errors + if (oldhelix == null) { + + if (_debug) + System.out.println(this.getClass().getSimpleName() + ": no old helix do the fit without MS scattering map"); + + // Reset the stereo hit positions to their nominal value + for (HelicalTrackHit hit : hitlist) + if (hit instanceof HelicalTrackCross) + ((HelicalTrackCross) hit).resetTrackDirection(); + // Do the fit + _status = _fitter.fit(hitlist); + SaveFit(); + + // Retrieve the helix parameters from this fit and save them in the seed + oldhelix = _fitter.getFit(); + seed.setHelix(oldhelix); + + if (_debug) + System.out.printf("%s: fit succeeded, will be used as seed, with chi2=%.3f and helix:\n%s \n", this.getClass().getSimpleName(), oldhelix.chisqtot(), oldhelix.toString()); + + // Calculate the multiple scattering angles for this helix + try { + seed.setScatterAngles(_scattering.FindScatters(oldhelix)); + } catch (Exception e) { + System.err.println(e); + if(_debug) + { + e.printStackTrace(); + } + return false; + } + + if(_debug) { + System.out.printf("%s: after calculating the MS map it has %d size:\n",this.getClass().getSimpleName(),seed.getMSMap().size()); + for(Map.Entry<HelicalTrackHit, MultipleScatter> ms : seed.getMSMap().entrySet()) { + System.out.printf("%s: Hit at layer %d and position %s has MS drpdhi=%f and dz=%f\n", + this.getClass().getSimpleName(),ms.getKey().Layer(),ms.getKey().getCorrectedPosition().toString(),ms.getValue().drphi(),ms.getValue().dz()); + } + } + + } + + if (_debug) + System.out.printf("%s: update the stereo hit positions with the old helix:\n%s \n", this.getClass().getSimpleName(), oldhelix.toString()); + + // Update the stereo hit positions and covariance matrices + CorrectStereoHits(hitlist, oldhelix); + + // Do a helix fit including MS errors + if(_debug) { + System.out.printf("%s: do the helix fit including MS map this time: \n",this.getClass().getSimpleName()); + for(Map.Entry<HelicalTrackHit, MultipleScatter> ms : seed.getMSMap().entrySet()) { + System.out.printf("%s: Hit at layer %d and position %s has MS drpdhi=%f and dz=%f\n", + this.getClass().getSimpleName(),ms.getKey().Layer(),ms.getKey().getCorrectedPosition().toString(),ms.getValue().drphi(),ms.getValue().dz()); + } + } + _status = _fitter.fit(hitlist, seed.getMSMap(), oldhelix); + SaveFit(); + + // Retrieve and save the new helix fit + _helix = _fitter.getFit(); + seed.setHelix(_helix); + if ((_trackCheck != null) && (!_trackCheck.checkSeed(seed))) + return false; + + // Set the non-holonomic constraint chi square...don't worry about this for now-- MG 10/29/2014 + // _constrain.setConstraintChisq(strategy, _helix, hitlist); + // If the total chi square is below the cut, we have a successful fit + boolean success = _helix.chisqtot() <= strategy.getMaxChisq(); + if (!success) + if (_diag != null) + _diag.fireFailedChisqCut(seed); + +// // If fit was successful, set the new multiple scattering angles + if (success) { + seed.setScatterAngles(_scattering.FindScatters(_helix)); + if (_debug) + System.out.printf("%s: this fit was successful, chi2=%f with resulting helix paramaters:\n%s\n", this.getClass().getSimpleName(), _helix.chisqtot(), _helix.toString()); // System.out.printf("%s: updated MS map before returning from fitCandidate():\n",this.getClass().getSimpleName()); + // for(Map.Entry<HelicalTrackHit, MultipleScatter> ms : seed.getMSMap().entrySet()) { + // System.out.printf("%s: Hit at layer %d and position %s has MS drpdhi=%f and dz=%f\n", + // this.getClass().getSimpleName(),ms.getKey().Layer(),ms.getKey().getCorrectedPosition().toString(),ms.getValue().drphi(),ms.getValue().dz()); + + } else if (_debug) + System.out.printf("%s: this fit with chi2=%f failed!\n", this.getClass().getSimpleName(), _helix.chisqtot()); + + return success; + } + + public void setDiagnostics(ISeedTrackerDiagnostics d) { + _diag = d; + return; + } + + public HelicalTrackFit getHelix() { + return _helix; + } + + public FitStatus getFitStatus() { + return _status; + } + + public CircleFit getCircleFit() { + return _circlefit; + } + + public SlopeInterceptLineFit getLineFit() { + return _linefit; + } + + public ZSegmentFit getZSegmentFit() { + return _zsegmentfit; + } + + public void setBField(double bfield) { + _bfield = bfield; + _scattering.setBField(_bfield); + // _constrain.setBField(_bfield); + return; + } + + public void setReferencePoint(double x, double y) { + // _fitter.setReferencePoint(x, y); + return; + } + + private void SaveFit() { + + // Default to no fit results when circle fit fails + _circlefit = null; + _linefit = null; + _zsegmentfit = null; + + // If we have a circle fit, try to save the line fit / zsegment fit results + if (_status == FitStatus.InconsistentSeed) + return; + if (_status == FitStatus.XYLineFitFailed) + return; + if (_status == FitStatus.SZLineFitFailed) + return; + // _linefit = _fitter.getLineFit(); +// _zsegmentfit = _fitter.getZSegmentFit(); + + return; + } + + private void CorrectStereoHits(List<HelicalTrackHit> hitlist, HelicalTrackFit helix) { + + // Get the PathMap - used to find the track direction at the hit location + Map<HelicalTrackHit, Double> pathmap = helix.PathMap(); + + // Loop over the hits and look for stereo hits + for (HelicalTrackHit hit : hitlist) + if (hit instanceof HelicalTrackCross) + if (pathmap.containsKey(hit)) { + // Found a stereo hit - calculate the track direction and pass it to the hit + ((HelicalTrackCross) hit).setTrackDirection(helix); + } + return; + } + + public void setDebug(boolean debug) { + this._debug = debug; + _scattering.setDebug(true); + } +} Copied: java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTrackReconDriver.java (from r1321, java/trunk/tracking/src/main/java/org/hps/recon/tracking/TrackerReconDriver.java) ============================================================================= --- java/trunk/tracking/src/main/java/org/hps/recon/tracking/TrackerReconDriver.java (original) +++ java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTrackReconDriver.java Thu Oct 30 08:16:13 2014 @@ -1,9 +1,9 @@ -package org.hps.recon.tracking; +package org.hps.recon.tracking.nobfield; import hep.physics.vec.BasicHep3Vector; - import java.util.List; - +import org.hps.recon.tracking.CoordinateTransformations; +import org.hps.recon.tracking.nobfield.StraightTracker; import org.lcsim.event.EventHeader; import org.lcsim.event.Track; import org.lcsim.event.base.BaseTrack; @@ -21,7 +21,7 @@ * * @author Matt Graham */ -public final class TrackerReconDriver extends Driver { +public final class StraightTrackReconDriver extends Driver { // Debug flag. private boolean debug = false; @@ -34,9 +34,9 @@ // Default B-field value. private double bfield = 0.5; // Tracking strategies resource path. - private String strategyResource = "HPS-Test-4pt1.xml"; + private String strategyResource = "HPS-Full.xml"; // Output track collection. - private String trackCollectionName = "MatchedTracks"; + private String trackCollectionName = "StraightTracks"; // HelicalTrackHit input collection. private String stInputCollectionName = "RotatedHelicalTrackHits"; // Include MS (done by removing XPlanes from the material manager results) @@ -49,7 +49,7 @@ private boolean _applySectorBinning = true; private double rmsTimeCut = -1; - public TrackerReconDriver() { + public StraightTrackReconDriver() { } public void setDebug(boolean debug) { @@ -134,15 +134,15 @@ */ private void initialize() { - // + // // 1) Driver to run Seed Tracker. // if (!strategyResource.startsWith("/")) { strategyResource = "/org/hps/recon/tracking/strategies/" + strategyResource; } List<SeedStrategy> sFinallist = StrategyXMLUtils.getStrategyListFromInputStream(this.getClass().getResourceAsStream(strategyResource)); - SeedTracker stFinal = new SeedTracker(sFinallist, this._useHPSMaterialManager, this.includeMS); - stFinal.setApplySectorBinning(_applySectorBinning); + StraightTracker stFinal = new StraightTracker(sFinallist, this._useHPSMaterialManager, this.includeMS); +// stFinal.setApplySectorBinning(_applySectorBinning); stFinal.setUseDefaultXPlane(false); stFinal.setDebug(this.debug); stFinal.setIterativeConfirmed(_iterativeConfirmed); @@ -156,10 +156,10 @@ // stFinal.setSectorParams(false); //this doesn't actually seem to do anything stFinal.setSectorParams(1, 10000); add(stFinal); - - if (rmsTimeCut > 0) { - stFinal.setTrackCheck(new HitTimeTrackCheck(rmsTimeCut)); - } +// +// if (rmsTimeCut > 0) { +// stFinal.setTrackCheck(new HitTimeTrackCheck(rmsTimeCut)); +// } } /** Added: java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTracker.java ============================================================================= --- java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTracker.java (added) +++ java/trunk/tracking/src/main/java/org/hps/recon/tracking/nobfield/StraightTracker.java Thu Oct 30 08:16:13 2014 @@ -0,0 +1,328 @@ +/* + * SeedTracker.java + * + * Created on August 16, 2005, 8:54 AM + * + */ +package org.hps.recon.tracking.nobfield; + +import hep.physics.vec.BasicHep3Vector; +import hep.physics.vec.Hep3Vector; + +import java.util.*; +import org.hps.recon.tracking.MaterialManager; +import org.hps.recon.tracking.MaterialSupervisor; +import org.lcsim.detector.ITransform3D; +import org.lcsim.event.EventHeader; +import org.lcsim.event.MCParticle; +import org.lcsim.fit.helicaltrack.HelicalTrackHit; +import org.lcsim.geometry.Detector; +import org.lcsim.recon.tracking.seedtracker.DefaultStrategy; +import org.lcsim.recon.tracking.seedtracker.HitManager; +import org.lcsim.recon.tracking.seedtracker.MakeTracks; +import org.lcsim.recon.tracking.seedtracker.SeedCandidate; +import org.lcsim.recon.tracking.seedtracker.SeedStrategy; +import org.lcsim.recon.tracking.seedtracker.TrackCheck; +import org.lcsim.recon.tracking.seedtracker.diagnostic.ISeedTrackerDiagnostics; +import org.lcsim.util.Driver; +import org.lcsim.util.aida.AIDA; + +/** + * Tracking algorithm based on forming track seeds from all 3-hit combinations, + * confirming this tentantive helix by requiring additional hits, and extending + * the track to additional layers. The operation of the algorithm is controlled + * by a list of SeedStrategy that define the tracker layers to be used and the + * cuts on the tracking algorithm. + * + * @author Richard Partridge + * @version 1.0 + */ +public class StraightTracker extends Driver { + + protected List<SeedStrategy> _strategylist; + protected ISeedTrackerDiagnostics _diag = null; + protected MaterialManager _materialmanager = new MaterialManager(); + protected HitManager _hitmanager; + protected StraightTrackFitter _helixfitter; + protected StraightTrackFinder _finder; + protected MakeTracks _maketracks; + protected Hep3Vector _IP = new BasicHep3Vector(0., 0., 0.); + protected double _bfield = 0.; + protected boolean _forceBField = false; + protected double _rtrk = 1000.; + protected boolean _autosectoring = false; + protected AIDA aida = AIDA.defaultInstance(); + protected boolean _timing = false; + protected String _inputCol = "HelicalTrackHits"; + private int _iterativeConfirmedFits = 0; + private boolean _debug = false; + + /** + * Creates a new instance of SeedTracker + */ + public StraightTracker() { + this(new DefaultStrategy().getStrategyList(), true, true); + } + + public StraightTracker(List<SeedStrategy> strategylist, boolean includeMS) { + initialize(strategylist, true, includeMS); + } + + public StraightTracker(List<SeedStrategy> strategylist, boolean useHPSMaterialManager, boolean includeMS) { + initialize(strategylist, useHPSMaterialManager, includeMS); + } + + private void initialize(List<SeedStrategy> strategylist, boolean useHPSMaterialManager, boolean includeMS) { + _strategylist = strategylist; + // Instantiate the hit manager + _hitmanager = new HitManager(); + // Instantiate the material manager for HPS, the helix fitter and seed track finder as tey depends on the material manager + if (useHPSMaterialManager) { + MaterialSupervisor materialSupervisor = new MaterialSupervisor(includeMS); + _materialmanager = materialSupervisor; + _helixfitter = new StraightTrackFitter(materialSupervisor); + } else { + MaterialManager materialmanager = new MaterialManager(includeMS); + _materialmanager = materialmanager; //mess around with types here... + _helixfitter = new StraightTrackFitter(materialmanager); + } + // Instantiate the Seed Finder + _finder = new StraightTrackFinder(_hitmanager, _helixfitter); + // Instantiate the Track Maker + _maketracks = new MakeTracks(); + } + + /** + * Set to enable debug output + * + * @param debug switch + */ + public void setDebug(boolean debug) { + _debug = true; + _materialmanager.setDebug(debug); + _helixfitter.setDebug(debug); + _finder.setDebug(debug); + } + + @Override + protected void process(EventHeader event) { + +// System.out.println("New event"); + // Pass the event to the diagnostics package + if (_diag != null) + _diag.setEvent(event); + + // Initialize timing + long last_time = System.currentTimeMillis(); + + // Get the hit collection from the event + List<HelicalTrackHit> hitcol = event.get(HelicalTrackHit.class, _inputCol); + if (_debug) + System.out.println("In " + this.getClass().getSimpleName() + ": Number of HelicalTrackHits=" + hitcol.size()); + // Sort the hits for this event + _hitmanager.OrganizeHits(hitcol); + + // Make the timing plots if requested + long start_time = System.currentTimeMillis(); + double dtime = ((double) (start_time - last_time)) / 1000.; + last_time = start_time; + if (_timing) + aida.cloud1D("Organize Hits").fill(dtime); + + // Make sure that we have cleared the list of track seeds in the finder + _finder.clearTrackSeedList(); + + // Loop over strategies and perform track finding + for (SeedStrategy strategy : _strategylist) { + + // Set the strategy for the diagnostics + if (_diag != null) + _diag.fireStrategyChanged(strategy); + + // Perform track finding under this strategy + _finder.FindTracks(strategy, _bfield); + if (_debug) + System.out.println("In " + this.getClass().getSimpleName() + ": Number of Tracks Found=" + _finder.getTrackSeeds().size()); + // Make the timing plots if requested + long time = System.currentTimeMillis(); + dtime = ((double) (time - last_time)) / 1000.; + last_time = time; + if (_timing) + aida.cloud1D("Tracking time for strategy " + strategy.getName()).fill(dtime); + } + + // Get the list of final list of SeedCandidates + List<SeedCandidate> trackseeds = _finder.getTrackSeeds(); + + if (_iterativeConfirmedFits > 0) { + // Iteratively re-fit tracks to take into account helix and hit position correlations + if (_debug) + System.out.printf("%s: Iteratively improve %d seeds\n", this.getClass().getSimpleName(), trackseeds.size()); + List<SeedCandidate> seedsToRemove = new ArrayList<SeedCandidate>(); + for (SeedCandidate seed : trackseeds) { + SeedStrategy strategy = seed.getSeedStrategy(); + boolean success = false; + for (int iterFit = 0; iterFit < _iterativeConfirmedFits; ++iterFit) + success = _helixfitter.FitCandidate(seed, strategy); + if (!success) + seedsToRemove.add(seed); +// else if (_debug) +// System.out.printf("%s: done iterating, this seed will be added to event:\n%s\n", this.getClass().getSimpleName(), seed.toString()); + } + for (SeedCandidate badseed : seedsToRemove) + trackseeds.remove(badseed); + } + + // Make tracks from the final list of track seeds + _maketracks.Process(event, trackseeds, _bfield); + + // Save the MC Particles that have been seeded / confirmed if diagnostics are enabled + if (_diag != null) { + Set<MCParticle> seededmcpset = _finder.getSeededMCParticles(); + List<MCParticle> seededmcp = new ArrayList<MCParticle>(seededmcpset); + event.put("SeededMCParticles", seededmcp, MCParticle.class, 0); + Set<MCParticle> confirmedmcpset = _finder.getConfirmedMCParticles(); + List<MCParticle> confirmedmcp = new ArrayList<MCParticle>(confirmedmcpset); + event.put("ConfirmedMCParticles", confirmedmcp, MCParticle.class, 0); + } + + // Clear the list of track seeds accumulated in the track finder + _finder.clearTrackSeedList(); + + // Make the total time plot if requested + long end_time = System.currentTimeMillis(); + dtime = ((double) (end_time - start_time)) / 1000.; + if (_timing) + aida.cloud1D("Total tracking time").fill(dtime); + + return; + } + + @Override + protected void detectorChanged(Detector detector) { + + // Only build the model when the detector is changed + _materialmanager.buildModel(detector); + + // Find the bfield and pass it to the helix fitter and diagnostic package + if (!_forceBField) + _bfield = detector.getFieldMap().getField(_IP).z(); + + if (_diag != null) + _diag.setBField(_bfield); + _helixfitter.setBField(_bfield); + + // Get the tracking radius + _rtrk = _materialmanager.getRMax(); + + // Set the sectoring parameters + if (_autosectoring) + _hitmanager.setSectorParams(_strategylist, _bfield, _rtrk); + } + + /** + * Specifiy the strategies to be used by the SeedTracker algorithm. Invoking + * this + * method will override the default strategies defined by the + * DefaultStrategy + * class. + * + * @param strategylist List of strategies to be used + */ + public void putStrategyList(List<SeedStrategy> strategylist) { + + // Save the strategy list + _strategylist = strategylist; + + // Set the sectoring parameters + if (_autosectoring) + _hitmanager.setSectorParams(strategylist, _bfield, _rtrk); + + return; + } + + public void setSectorParams(int nphi, double dz) { + _hitmanager.setSectorParams(nphi, dz); + _autosectoring = false; + return; + } + + public void setDiagnostics(ISeedTrackerDiagnostics d) { + + // Set the diagnostic package + _diag = d; + _helixfitter.setDiagnostics(_diag); + _finder.setDiagnostic(_diag); + + // Pass the hit manager, material manager, and bfield to the diagnostic package + if (_diag != null) { + _diag.setMaterialManager(_materialmanager); + _diag.setHitManager(_hitmanager); + _diag.setBField(_bfield); + } + } + + public void setTimingPlots(boolean timing) { + _timing = timing; + } + + public void setTrkCollectionName(String name) { + _maketracks.setTrkCollectionName(name); + } + + public void setInputCollectionName(String name) { + _inputCol = name; + } + + public void setMaterialManagerTransform(ITransform3D _detToTrk) { + _materialmanager.setTransform(_detToTrk); + } + + /** + * Set the maximum number of fits used to confirm or extend a seed. + * + * @param maxfit maximum number of fits + */ + public void setMaxFit(int maxfit) { + _finder.setMaxFit(maxfit); + } + + public void setBField(double bfield) { + _forceBField = true; + _bfield = bfield; + } + + public void setReferencePoint(double xref, double yref) { + _helixfitter.setReferencePoint(xref, yref); + } + + public void setSectorParams(boolean sector) { + _hitmanager.setDoSectoring(sector); + } + + /** + * Set {@link TrackCheck} object to be used by the track finding algorithm. + * If this method is never called, no external checking of seeds and tracks + * is performed. + */ + public void setTrackCheck(TrackCheck trackCheck) { + _finder._trackCheck = trackCheck; + _maketracks._trackCheck = trackCheck; + } + + /** + * Set the maximum number of iterative fits on a confirmed/extended + * candidate. + * + * @param maxfits maximum number of fits + */ + public void setIterativeConfirmed(int maxfits) { + this._iterativeConfirmedFits = maxfits; + } + + public void setUseDefaultXPlane(boolean useDefault) { + _materialmanager.setDefaultXPlaneUsage(useDefault); + + } + +}