Examples with ImageSource gdsc.smlm.results.ImageSource used on opensource projects

Example 1 with ImageSource

use of gdsc.smlm.results.ImageSource in project GDSC-SMLM by aherbert.

the class PSFEstimator method calculateStatistics.

private boolean calculateStatistics(PeakFit fitter, double[] params, double[] params_dev) {
    debug("  Fitting PSF");
    swapStatistics();
    // Create the fit engine using the PeakFit plugin
    FitConfiguration fitConfig = config.getFitConfiguration();
    fitConfig.setInitialAngle((float) params[0]);
    fitConfig.setInitialPeakStdDev0((float) params[1]);
    fitConfig.setInitialPeakStdDev1((float) params[2]);
    ImageStack stack = imp.getImageStack();
    Rectangle roi = stack.getProcessor(1).getRoi();
    ImageSource source = new IJImageSource(imp);
    // Allow interlaced data by wrapping the image source
    if (interlacedData) {
        source = new InterlacedImageSource(source, dataStart, dataBlock, dataSkip);
    }
    // Allow frame aggregation by wrapping the image source
    if (integrateFrames > 1) {
        source = new AggregatedImageSource(source, integrateFrames);
    }
    fitter.initialiseImage(source, roi, true);
    fitter.addPeakResults(this);
    fitter.initialiseFitting();
    FitEngine engine = fitter.createFitEngine();
    // Use random slices
    int[] slices = new int[stack.getSize()];
    for (int i = 0; i < slices.length; i++) slices[i] = i + 1;
    Random rand = new Random();
    rand.shuffle(slices);
    IJ.showStatus("Fitting ...");
    // Use multi-threaded code for speed
    int i;
    for (i = 0; i < slices.length; i++) {
        int slice = slices[i];
        //debug("  Processing slice = %d\n", slice);
        IJ.showProgress(size(), settings.numberOfPeaks);
        ImageProcessor ip = stack.getProcessor(slice);
        // stack processor does not set the bounds required by ImageConverter
        ip.setRoi(roi);
        FitJob job = new FitJob(slice, ImageConverter.getData(ip), roi);
        engine.run(job);
        if (sampleSizeReached() || Utils.isInterrupted()) {
            break;
        }
    }
    if (Utils.isInterrupted()) {
        IJ.showProgress(1);
        engine.end(true);
        return false;
    }
    // Wait until we have enough results
    while (!sampleSizeReached() && !engine.isQueueEmpty()) {
        IJ.showProgress(size(), settings.numberOfPeaks);
        try {
            Thread.sleep(50);
        } catch (InterruptedException e) {
            break;
        }
    }
    // End now if we have enough samples
    engine.end(sampleSizeReached());
    IJ.showStatus("");
    IJ.showProgress(1);
    // This count will be an over-estimate given that the provider is ahead of the consumer
    // in this multi-threaded system
    debug("  Processed %d/%d slices (%d peaks)", i, slices.length, size());
    setParams(ANGLE, params, params_dev, sampleNew[ANGLE]);
    setParams(X, params, params_dev, sampleNew[X]);
    setParams(Y, params, params_dev, sampleNew[Y]);
    if (settings.showHistograms) {
        int[] idList = new int[NAMES.length];
        int count = 0;
        boolean requireRetile = false;
        for (int ii = 0; ii < 3; ii++) {
            if (sampleNew[ii].getN() == 0)
                continue;
            StoredDataStatistics stats = new StoredDataStatistics(sampleNew[ii].getValues());
            idList[count++] = Utils.showHistogram(TITLE, stats, NAMES[ii], 0, 0, settings.histogramBins, "Mean = " + Utils.rounded(stats.getMean()) + ". Median = " + Utils.rounded(sampleNew[ii].getPercentile(50)));
            requireRetile = requireRetile || Utils.isNewWindow();
        }
        if (requireRetile && count > 0) {
            new WindowOrganiser().tileWindows(Arrays.copyOf(idList, count));
        }
    }
    if (size() < 2) {
        log("ERROR: Insufficient number of fitted peaks, terminating ...");
        return false;
    }
    return true;
}

Example 2 with ImageSource

use of gdsc.smlm.results.ImageSource in project GDSC-SMLM by aherbert.

the class PeakFit method setup.

/*
	 * (non-Javadoc)
	 * 
	 * @see ij.plugin.filter.PlugInFilter#setup(java.lang.String, ij.ImagePlus)
	 */
public int setup(String arg, ImagePlus imp) {
    SMLMUsageTracker.recordPlugin(this.getClass(), arg);
    plugin_flags = FLAGS;
    extraOptions = Utils.isExtraOptions();
    maximaIdentification = (arg != null && arg.contains("spot"));
    fitMaxima = (arg != null && arg.contains("maxima"));
    simpleFit = (arg != null && arg.contains("simple"));
    boolean runSeries = (arg != null && arg.contains("series"));
    ImageSource imageSource = null;
    if (fitMaxima) {
        imp = null;
        // The image source will be found from the peak results.
        if (!showMaximaDialog())
            return DONE;
        MemoryPeakResults results = ResultsManager.loadInputResults(inputOption, false);
        if (results == null || results.size() == 0) {
            IJ.error(TITLE, "No results could be loaded");
            return DONE;
        }
        // Check for single frame
        singleFrame = results.getHead().getFrame();
        for (PeakResult result : results.getResults()) {
            if (singleFrame != result.getFrame()) {
                singleFrame = 0;
                break;
            }
        }
        imageSource = results.getSource();
        plugin_flags |= NO_IMAGE_REQUIRED;
    } else if (runSeries) {
        imp = null;
        // Select input folder
        String inputDirectory;
        inputDirectory = IJ.getDirectory("Select image series ...");
        //inputDirectory = getInputDirectory("Select image series ...");
        if (inputDirectory == null)
            return DONE;
        // Load input series ...
        SeriesOpener series;
        if (extraOptions)
            series = new SeriesOpener(inputDirectory, true, numberOfThreads);
        else
            series = new SeriesOpener(inputDirectory);
        if (series.getNumberOfImages() == 0) {
            IJ.error(TITLE, "No images in the selected directory:\n" + inputDirectory);
            return DONE;
        }
        SeriesImageSource seriesImageSource = new SeriesImageSource(getName(series.getImageList()), series);
        seriesImageSource.setLogProgress(true);
        if (extraOptions) {
            numberOfThreads = Math.max(1, series.getNumberOfThreads());
            seriesImageSource.setNumberOfThreads(numberOfThreads);
        }
        imageSource = seriesImageSource;
        plugin_flags |= NO_IMAGE_REQUIRED;
    } else {
        if (imp == null) {
            IJ.noImage();
            return DONE;
        }
        // Check it is not a previous result
        if (imp.getTitle().endsWith(IJImagePeakResults.IMAGE_SUFFIX)) {
            IJImageSource tmpImageSource = null;
            // Check the image to see if it has an image source XML structure in the info property
            Object o = imp.getProperty("Info");
            Pattern pattern = Pattern.compile("Source: (<.*IJImageSource>.*<.*IJImageSource>)", Pattern.DOTALL);
            Matcher match = pattern.matcher((o == null) ? "" : o.toString());
            if (match.find()) {
                ImageSource source = ImageSource.fromXML(match.group(1));
                if (source instanceof IJImageSource) {
                    tmpImageSource = (IJImageSource) source;
                    if (!tmpImageSource.open()) {
                        tmpImageSource = null;
                    } else {
                        imp = WindowManager.getImage(tmpImageSource.getName());
                    }
                }
            }
            if (tmpImageSource == null) {
                // Look for a parent using the title
                String parentTitle = imp.getTitle().substring(0, imp.getTitle().length() - IJImagePeakResults.IMAGE_SUFFIX.length() - 1);
                ImagePlus parentImp = WindowManager.getImage(parentTitle);
                if (parentImp != null) {
                    tmpImageSource = new IJImageSource(parentImp);
                    imp = parentImp;
                }
            }
            String message = "The selected image may be a previous fit result";
            if (tmpImageSource != null) {
                // are missing
                if (!Utils.isNullOrEmpty(tmpImageSource.getName()))
                    message += " of: \n \n" + tmpImageSource.getName();
                message += " \n \nFit the parent?";
            } else
                message += " \n \nDo you want to continue?";
            YesNoCancelDialog d = new YesNoCancelDialog(null, TITLE, message);
            if (tmpImageSource == null) {
                if (!d.yesPressed())
                    return DONE;
            } else {
                if (d.yesPressed())
                    imageSource = tmpImageSource;
                if (d.cancelPressed())
                    return DONE;
            }
        }
        if (imageSource == null)
            imageSource = new IJImageSource(imp);
    }
    time = -1;
    if (!initialiseImage(imageSource, getBounds(imp), false)) {
        IJ.error(TITLE, "Failed to initialise the source image: " + imageSource.getName());
        return DONE;
    }
    int flags = showDialog(imp);
    if ((flags & DONE) == 0) {
        // Repeat so that we pass in the selected option for ignoring the bounds.
        // This should not be necessary since it is set within the readDialog method.
        //if (ignoreBoundsForNoise)
        //	initialiseImage(imageSource, bounds, ignoreBoundsForNoise);
        initialiseFitting();
    }
    return flags;
}

Example 3 with ImageSource

use of gdsc.smlm.results.ImageSource in project GDSC-SMLM by aherbert.

the class TraceMolecules method fitTraces.

private void fitTraces(MemoryPeakResults results, Trace[] traces) {
    // Check if the original image is open and the fit configuration can be extracted
    ImageSource source = results.getSource();
    if (source == null)
        return;
    if (!source.open())
        return;
    FitEngineConfiguration config = (FitEngineConfiguration) XmlUtils.fromXML(results.getConfiguration());
    if (config == null)
        return;
    // Show a dialog asking if the traces should be refit
    ExtendedGenericDialog gd = new ExtendedGenericDialog(TITLE);
    gd.addMessage("Do you want to fit the traces as a single peak using a combined image?");
    gd.addCheckbox("Fit_closest_to_centroid", !fitOnlyCentroid);
    gd.addSlider("Distance_threshold", 0.01, 3, distanceThreshold);
    gd.addSlider("Expansion_factor", 1, 4.5, expansionFactor);
    // Allow fitting settings to be adjusted
    FitConfiguration fitConfig = config.getFitConfiguration();
    gd.addMessage("--- Gaussian fitting ---");
    String[] filterTypes = SettingsManager.getNames((Object[]) DataFilterType.values());
    gd.addChoice("Spot_filter_type", filterTypes, filterTypes[config.getDataFilterType().ordinal()]);
    String[] filterNames = SettingsManager.getNames((Object[]) DataFilter.values());
    gd.addChoice("Spot_filter", filterNames, filterNames[config.getDataFilter(0).ordinal()]);
    gd.addSlider("Smoothing", 0, 2.5, config.getSmooth(0));
    gd.addSlider("Search_width", 0.5, 2.5, config.getSearch());
    gd.addSlider("Border", 0.5, 2.5, config.getBorder());
    gd.addSlider("Fitting_width", 2, 4.5, config.getFitting());
    String[] solverNames = SettingsManager.getNames((Object[]) FitSolver.values());
    gd.addChoice("Fit_solver", solverNames, solverNames[fitConfig.getFitSolver().ordinal()]);
    String[] functionNames = SettingsManager.getNames((Object[]) FitFunction.values());
    gd.addChoice("Fit_function", functionNames, functionNames[fitConfig.getFitFunction().ordinal()]);
    String[] criteriaNames = SettingsManager.getNames((Object[]) FitCriteria.values());
    gd.addChoice("Fit_criteria", criteriaNames, criteriaNames[fitConfig.getFitCriteria().ordinal()]);
    gd.addNumericField("Significant_digits", fitConfig.getSignificantDigits(), 0);
    gd.addNumericField("Coord_delta", fitConfig.getDelta(), 4);
    gd.addNumericField("Lambda", fitConfig.getLambda(), 4);
    gd.addNumericField("Max_iterations", fitConfig.getMaxIterations(), 0);
    gd.addNumericField("Fail_limit", config.getFailuresLimit(), 0);
    gd.addCheckbox("Include_neighbours", config.isIncludeNeighbours());
    gd.addSlider("Neighbour_height", 0.01, 1, config.getNeighbourHeightThreshold());
    gd.addSlider("Residuals_threshold", 0.01, 1, config.getResidualsThreshold());
    //gd.addSlider("Duplicate_distance", 0, 1.5, fitConfig.getDuplicateDistance());
    gd.addMessage("--- Peak filtering ---\nDiscard fits that shift; are too low; or expand/contract");
    gd.addCheckbox("Smart_filter", fitConfig.isSmartFilter());
    gd.addCheckbox("Disable_simple_filter", fitConfig.isDisableSimpleFilter());
    gd.addSlider("Shift_factor", 0.01, 2, fitConfig.getCoordinateShiftFactor());
    gd.addNumericField("Signal_strength", fitConfig.getSignalStrength(), 2);
    gd.addNumericField("Min_photons", fitConfig.getMinPhotons(), 0);
    gd.addSlider("Min_width_factor", 0, 0.99, fitConfig.getMinWidthFactor());
    gd.addSlider("Width_factor", 1.01, 5, fitConfig.getWidthFactor());
    gd.addNumericField("Precision", fitConfig.getPrecisionThreshold(), 2);
    gd.addCheckbox("Debug_failures", debugFailures);
    gd.showDialog();
    if (!gd.wasOKed()) {
        source.close();
        return;
    }
    // Get parameters for the fit
    fitOnlyCentroid = !gd.getNextBoolean();
    distanceThreshold = (float) gd.getNextNumber();
    expansionFactor = (float) gd.getNextNumber();
    config.setDataFilterType(gd.getNextChoiceIndex());
    config.setDataFilter(gd.getNextChoiceIndex(), Math.abs(gd.getNextNumber()), 0);
    config.setSearch(gd.getNextNumber());
    config.setBorder(gd.getNextNumber());
    config.setFitting(gd.getNextNumber());
    fitConfig.setFitSolver(gd.getNextChoiceIndex());
    fitConfig.setFitFunction(gd.getNextChoiceIndex());
    fitConfig.setFitCriteria(gd.getNextChoiceIndex());
    fitConfig.setSignificantDigits((int) gd.getNextNumber());
    fitConfig.setDelta(gd.getNextNumber());
    fitConfig.setLambda(gd.getNextNumber());
    fitConfig.setMaxIterations((int) gd.getNextNumber());
    config.setFailuresLimit((int) gd.getNextNumber());
    config.setIncludeNeighbours(gd.getNextBoolean());
    config.setNeighbourHeightThreshold(gd.getNextNumber());
    config.setResidualsThreshold(gd.getNextNumber());
    fitConfig.setSmartFilter(gd.getNextBoolean());
    fitConfig.setDisableSimpleFilter(gd.getNextBoolean());
    fitConfig.setCoordinateShiftFactor(gd.getNextNumber());
    fitConfig.setSignalStrength(gd.getNextNumber());
    fitConfig.setMinPhotons(gd.getNextNumber());
    fitConfig.setMinWidthFactor(gd.getNextNumber());
    fitConfig.setWidthFactor(gd.getNextNumber());
    fitConfig.setPrecisionThreshold(gd.getNextNumber());
    // Check arguments
    try {
        Parameters.isAboveZero("Distance threshold", distanceThreshold);
        Parameters.isAbove("Expansion factor", expansionFactor, 1);
        Parameters.isAboveZero("Search_width", config.getSearch());
        Parameters.isAboveZero("Fitting_width", config.getFitting());
        Parameters.isAboveZero("Significant digits", fitConfig.getSignificantDigits());
        Parameters.isAboveZero("Delta", fitConfig.getDelta());
        Parameters.isAboveZero("Lambda", fitConfig.getLambda());
        Parameters.isAboveZero("Max iterations", fitConfig.getMaxIterations());
        Parameters.isPositive("Failures limit", config.getFailuresLimit());
        Parameters.isPositive("Neighbour height threshold", config.getNeighbourHeightThreshold());
        Parameters.isPositive("Residuals threshold", config.getResidualsThreshold());
        Parameters.isPositive("Coordinate Shift factor", fitConfig.getCoordinateShiftFactor());
        Parameters.isPositive("Signal strength", fitConfig.getSignalStrength());
        Parameters.isPositive("Min photons", fitConfig.getMinPhotons());
        Parameters.isPositive("Min width factor", fitConfig.getMinWidthFactor());
        Parameters.isPositive("Width factor", fitConfig.getWidthFactor());
        Parameters.isPositive("Precision threshold", fitConfig.getPrecisionThreshold());
    } catch (IllegalArgumentException e) {
        IJ.error(TITLE, e.getMessage());
        source.close();
        return;
    }
    debugFailures = gd.getNextBoolean();
    if (!PeakFit.configureSmartFilter(globalSettings, filename))
        return;
    if (!PeakFit.configureDataFilter(globalSettings, filename, false))
        return;
    if (!PeakFit.configureFitSolver(globalSettings, filename, false))
        return;
    // Adjust settings for a single maxima
    config.setIncludeNeighbours(false);
    fitConfig.setDuplicateDistance(0);
    // Create a fit engine
    MemoryPeakResults refitResults = new MemoryPeakResults();
    refitResults.copySettings(results);
    refitResults.setName(results.getName() + " Trace Fit");
    refitResults.setSortAfterEnd(true);
    refitResults.begin();
    // No border since we know where the peaks are and we must not miss them due to truncated searching 
    FitEngine engine = new FitEngine(config, refitResults, Prefs.getThreads(), FitQueue.BLOCKING);
    // Either : Only fit the centroid
    // or     : Extract a bigger region, allowing all fits to run as normal and then 
    //          find the correct spot using Euclidian distance.
    // Set up the limits
    final double stdDev = FastMath.max(fitConfig.getInitialPeakStdDev0(), fitConfig.getInitialPeakStdDev1());
    float fitWidth = (float) (stdDev * config.getFitting() * ((fitOnlyCentroid) ? 1 : expansionFactor));
    IJ.showStatus("Refitting traces ...");
    List<JobItem> jobItems = new ArrayList<JobItem>(traces.length);
    int singles = 0;
    int fitted = 0;
    for (int n = 0; n < traces.length; n++) {
        Trace trace = traces[n];
        if (n % 32 == 0)
            IJ.showProgress(n, traces.length);
        // Skip traces with one peak
        if (trace.size() == 1) {
            singles++;
            // Use the synchronized method to avoid thread clashes with the FitEngine
            refitResults.addSync(trace.getHead());
            continue;
        }
        Rectangle bounds = new Rectangle();
        double[] combinedNoise = new double[1];
        float[] data = buildCombinedImage(source, trace, fitWidth, bounds, combinedNoise, false);
        if (data == null)
            continue;
        // Fit the combined image
        FitParameters params = new FitParameters();
        params.noise = (float) combinedNoise[0];
        float[] centre = trace.getCentroid();
        if (fitOnlyCentroid) {
            int newX = (int) Math.round(centre[0]) - bounds.x;
            int newY = (int) Math.round(centre[1]) - bounds.y;
            params.maxIndices = new int[] { newY * bounds.width + newX };
        } else {
            params.filter = new ArrayList<float[]>();
            params.filter.add(new float[] { centre[0] - bounds.x, centre[1] - bounds.y });
            params.distanceThreshold = distanceThreshold;
        }
        // This is not needed since the bounds are passed using the FitJob
        //params.setOffset(new float[] { bounds.x, bounds.y });
        int startT = trace.getHead().getFrame();
        params.endT = trace.getTail().getFrame();
        ParameterisedFitJob job = new ParameterisedFitJob(n, params, startT, data, bounds);
        jobItems.add(new JobItem(job, trace, centre));
        engine.run(job);
        fitted++;
    }
    engine.end(false);
    IJ.showStatus("");
    IJ.showProgress(1);
    // Check the success ...
    FitStatus[] values = FitStatus.values();
    int[] statusCount = new int[values.length + 1];
    ArrayList<String> names = new ArrayList<String>(Arrays.asList(SettingsManager.getNames((Object[]) values)));
    names.add(String.format("No maxima within %.2f of centroid", distanceThreshold));
    int separated = 0;
    int success = 0;
    final int debugLimit = 3;
    for (JobItem jobItem : jobItems) {
        int id = jobItem.getId();
        ParameterisedFitJob job = jobItem.job;
        Trace trace = jobItem.trace;
        int[] indices = job.getIndices();
        FitResult fitResult = null;
        int status;
        if (indices.length < 1) {
            status = values.length;
        } else if (indices.length > 1) {
            // Choose the first OK result. This is all that matters for the success reporting
            for (int n = 0; n < indices.length; n++) {
                if (job.getFitResult(n).getStatus() == FitStatus.OK) {
                    fitResult = job.getFitResult(n);
                    break;
                }
            }
            // Otherwise use the closest failure. 
            if (fitResult == null) {
                final float[] centre = traces[id].getCentroid();
                double minD = Double.POSITIVE_INFINITY;
                for (int n = 0; n < indices.length; n++) {
                    // Since the fit has failed we use the initial parameters
                    final double[] params = job.getFitResult(n).getInitialParameters();
                    final double dx = params[Gaussian2DFunction.X_POSITION] - centre[0];
                    final double dy = params[Gaussian2DFunction.Y_POSITION] - centre[1];
                    final double d = dx * dx + dy * dy;
                    if (minD > d) {
                        minD = d;
                        fitResult = job.getFitResult(n);
                    }
                }
            }
            status = fitResult.getStatus().ordinal();
        } else {
            fitResult = job.getFitResult(0);
            status = fitResult.getStatus().ordinal();
        }
        // All jobs have only one peak
        statusCount[status]++;
        // Debug why any fits failed
        if (fitResult == null || fitResult.getStatus() != FitStatus.OK) {
            refitResults.addAll(trace.getPoints());
            separated += trace.size();
            if (debugFailures) {
                FitStatus s = (fitResult == null) ? FitStatus.UNKNOWN : fitResult.getStatus();
                // Only display the first n per category to limit the number of images
                double[] noise = new double[1];
                if (statusCount[status] <= debugLimit) {
                    Rectangle bounds = new Rectangle();
                    buildCombinedImage(source, trace, fitWidth, bounds, noise, true);
                    float[] centre = trace.getCentroid();
                    Utils.display(String.format("Trace %d (n=%d) : x=%f,y=%f", id, trace.size(), centre[0], centre[1]), slices);
                    switch(s) {
                        case INSUFFICIENT_PRECISION:
                            float precision = (Float) fitResult.getStatusData();
                            IJ.log(String.format("Trace %d (n=%d) : %s = %f", id, trace.size(), names.get(status), precision));
                            break;
                        case INSUFFICIENT_SIGNAL:
                            if (noise[0] == 0)
                                noise[0] = getCombinedNoise(trace);
                            float snr = (Float) fitResult.getStatusData();
                            IJ.log(String.format("Trace %d (n=%d) : %s = %f (noise=%.2f)", id, trace.size(), names.get(status), snr, noise[0]));
                            break;
                        case COORDINATES_MOVED:
                        case OUTSIDE_FIT_REGION:
                        case WIDTH_DIVERGED:
                            float[] shift = (float[]) fitResult.getStatusData();
                            IJ.log(String.format("Trace %d (n=%d) : %s = %.3f,%.3f", id, trace.size(), names.get(status), shift[0], shift[1]));
                            break;
                        default:
                            IJ.log(String.format("Trace %d (n=%d) : %s", id, trace.size(), names.get(status)));
                            break;
                    }
                }
            }
        } else {
            success++;
            if (debugFailures) {
                // Only display the first n per category to limit the number of images
                double[] noise = new double[1];
                if (statusCount[status] <= debugLimit) {
                    Rectangle bounds = new Rectangle();
                    buildCombinedImage(source, trace, fitWidth, bounds, noise, true);
                    float[] centre = trace.getCentroid();
                    Utils.display(String.format("Trace %d (n=%d) : x=%f,y=%f", id, trace.size(), centre[0], centre[1]), slices);
                }
            }
        }
    }
    IJ.log(String.format("Trace fitting : %d singles : %d / %d fitted : %d separated", singles, success, fitted, separated));
    if (separated > 0) {
        IJ.log("Reasons for fit failure :");
        // Start at i=1 to skip FitStatus.OK
        for (int i = 1; i < statusCount.length; i++) {
            if (statusCount[i] != 0)
                IJ.log("  " + names.get(i) + " = " + statusCount[i]);
        }
    }
    refitResults.end();
    MemoryPeakResults.addResults(refitResults);
    source.close();
}

Example 4 with ImageSource

use of gdsc.smlm.results.ImageSource in project GDSC-SMLM by aherbert.

the class SpotInspector method run.

/*
	 * (non-Javadoc)
	 * 
	 * @see ij.plugin.PlugIn#run(java.lang.String)
	 */
public void run(String arg) {
    SMLMUsageTracker.recordPlugin(this.getClass(), arg);
    if (MemoryPeakResults.isMemoryEmpty()) {
        IJ.error(TITLE, "No localisations in memory");
        return;
    }
    if (!showDialog())
        return;
    // Load the results
    results = ResultsManager.loadInputResults(inputOption, false);
    if (results == null || results.size() == 0) {
        IJ.error(TITLE, "No results could be loaded");
        IJ.showStatus("");
        return;
    }
    // Check if the original image is open
    ImageSource source = results.getSource();
    if (source == null) {
        IJ.error(TITLE, "Unknown original source image");
        return;
    }
    source = source.getOriginal();
    if (!source.open()) {
        IJ.error(TITLE, "Cannot open original source image: " + source.toString());
        return;
    }
    final float stdDevMax = getStandardDeviation(results);
    if (stdDevMax < 0) {
        // TODO - Add dialog to get the initial peak width
        IJ.error(TITLE, "Fitting configuration (for initial peak width) is not available");
        return;
    }
    // Rank spots
    rankedResults = new ArrayList<PeakResultRank>(results.size());
    final double a = results.getNmPerPixel();
    final double gain = results.getGain();
    final boolean emCCD = results.isEMCCD();
    for (PeakResult r : results.getResults()) {
        float[] score = getScore(r, a, gain, emCCD, stdDevMax);
        rankedResults.add(new PeakResultRank(r, score[0], score[1]));
    }
    Collections.sort(rankedResults);
    // Prepare results table. Get bias if necessary
    if (showCalibratedValues) {
        // Get a bias if required
        Calibration calibration = results.getCalibration();
        if (calibration.getBias() == 0) {
            ExtendedGenericDialog gd = new ExtendedGenericDialog(TITLE);
            gd.addMessage("Calibrated results requires a camera bias");
            gd.addNumericField("Camera_bias (ADUs)", calibration.getBias(), 2);
            gd.showDialog();
            if (!gd.wasCanceled()) {
                calibration.setBias(Math.abs(gd.getNextNumber()));
            }
        }
    }
    IJTablePeakResults table = new IJTablePeakResults(false, results.getName(), true);
    table.copySettings(results);
    table.setTableTitle(TITLE);
    table.setAddCounter(true);
    table.setShowCalibratedValues(showCalibratedValues);
    table.begin();
    // Add a mouse listener to jump to the frame for the clicked line
    textPanel = table.getResultsWindow().getTextPanel();
    // We must ignore old instances of this class from the mouse listeners
    id = ++currentId;
    textPanel.addMouseListener(this);
    // Add results to the table
    int n = 0;
    for (PeakResultRank rank : rankedResults) {
        rank.rank = n++;
        PeakResult r = rank.peakResult;
        table.add(r.getFrame(), r.origX, r.origY, r.origValue, r.error, r.noise, r.params, r.paramsStdDev);
    }
    table.end();
    if (plotScore || plotHistogram) {
        // Get values for the plots
        float[] xValues = null, yValues = null;
        double yMin, yMax;
        int spotNumber = 0;
        xValues = new float[rankedResults.size()];
        yValues = new float[xValues.length];
        for (PeakResultRank rank : rankedResults) {
            xValues[spotNumber] = spotNumber + 1;
            yValues[spotNumber++] = recoverScore(rank.score);
        }
        // Set the min and max y-values using 1.5 x IQR 
        DescriptiveStatistics stats = new DescriptiveStatistics();
        for (float v : yValues) stats.addValue(v);
        if (removeOutliers) {
            double lower = stats.getPercentile(25);
            double upper = stats.getPercentile(75);
            double iqr = upper - lower;
            yMin = FastMath.max(lower - iqr, stats.getMin());
            yMax = FastMath.min(upper + iqr, stats.getMax());
            IJ.log(String.format("Data range: %f - %f. Plotting 1.5x IQR: %f - %f", stats.getMin(), stats.getMax(), yMin, yMax));
        } else {
            yMin = stats.getMin();
            yMax = stats.getMax();
            IJ.log(String.format("Data range: %f - %f", yMin, yMax));
        }
        plotScore(xValues, yValues, yMin, yMax);
        plotHistogram(yValues, yMin, yMax);
    }
    // Extract spots into a stack
    final int w = source.getWidth();
    final int h = source.getHeight();
    final int size = 2 * radius + 1;
    ImageStack spots = new ImageStack(size, size, rankedResults.size());
    // To assist the extraction of data from the image source, process them in time order to allow 
    // frame caching. Then set the appropriate slice in the result stack
    Collections.sort(rankedResults, new Comparator<PeakResultRank>() {

        public int compare(PeakResultRank o1, PeakResultRank o2) {
            if (o1.peakResult.getFrame() < o2.peakResult.getFrame())
                return -1;
            if (o1.peakResult.getFrame() > o2.peakResult.getFrame())
                return 1;
            return 0;
        }
    });
    for (PeakResultRank rank : rankedResults) {
        PeakResult r = rank.peakResult;
        // Extract image
        // Note that the coordinates are relative to the middle of the pixel (0.5 offset)
        // so do not round but simply convert to int
        final int x = (int) (r.params[Gaussian2DFunction.X_POSITION]);
        final int y = (int) (r.params[Gaussian2DFunction.Y_POSITION]);
        // Extract a region but crop to the image bounds
        int minX = x - radius;
        int minY = y - radius;
        int maxX = FastMath.min(x + radius + 1, w);
        int maxY = FastMath.min(y + radius + 1, h);
        int padX = 0, padY = 0;
        if (minX < 0) {
            padX = -minX;
            minX = 0;
        }
        if (minY < 0) {
            padY = -minY;
            minY = 0;
        }
        int sizeX = maxX - minX;
        int sizeY = maxY - minY;
        float[] data = source.get(r.getFrame(), new Rectangle(minX, minY, sizeX, sizeY));
        // Prevent errors with missing data
        if (data == null)
            data = new float[sizeX * sizeY];
        ImageProcessor spotIp = new FloatProcessor(sizeX, sizeY, data, null);
        // Pad if necessary, i.e. the crop is too small for the stack
        if (padX > 0 || padY > 0 || sizeX < size || sizeY < size) {
            ImageProcessor spotIp2 = spotIp.createProcessor(size, size);
            spotIp2.insert(spotIp, padX, padY);
            spotIp = spotIp2;
        }
        int slice = rank.rank + 1;
        spots.setPixels(spotIp.getPixels(), slice);
        spots.setSliceLabel(Utils.rounded(rank.originalScore), slice);
    }
    source.close();
    ImagePlus imp = Utils.display(TITLE, spots);
    imp.setRoi((PointRoi) null);
    // Make bigger		
    for (int i = 10; i-- > 0; ) imp.getWindow().getCanvas().zoomIn(imp.getWidth() / 2, imp.getHeight() / 2);
}