Search in sources :

Example 1 with IJTablePeakResults

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

the class GaussianFit method runFinal.

/**
	 * Perform fitting using the chosen maxima. Update the overlay if successful.
	 * 
	 * @param ip
	 *            The input image
	 */
private void runFinal(ImageProcessor ip) {
    ip.reset();
    Rectangle bounds = ip.getRoi();
    // Crop to the ROI
    float[] data = ImageConverter.getData(ip);
    int width = bounds.width;
    int height = bounds.height;
    // Sort the maxima
    float[] smoothData = data;
    if (getSmooth() > 0) {
        // Smoothing destructively modifies the data so create a copy
        smoothData = Arrays.copyOf(data, width * height);
        AverageFilter filter = new AverageFilter();
        //filter.blockAverage(smoothData, width, height, smooth);
        if (smooth <= border)
            filter.stripedBlockAverageInternal(smoothData, width, height, (float) smooth);
        else
            filter.stripedBlockAverage(smoothData, width, height, (float) smooth);
    }
    Sort.sort(maxIndices, smoothData);
    // Show the candidate peaks
    if (maxIndices.length > 0) {
        String message = String.format("Identified %d peaks", maxIndices.length);
        if (isLogProgress()) {
            IJ.log(message);
            for (int index : maxIndices) {
                IJ.log(String.format("  %.2f @ [%d,%d]", data[index], bounds.x + index % width, bounds.y + index / width));
            }
        }
        // Check whether to run if the number of peaks is large
        if (maxIndices.length > 10) {
            GenericDialog gd = new GenericDialog("Warning");
            gd.addMessage(message + "\nDo you want to fit?");
            gd.showDialog();
            if (gd.wasCanceled())
                return;
        }
    } else {
        IJ.log("No maxima identified");
        return;
    }
    results = new IJTablePeakResults(showDeviations, imp.getTitle() + " [" + imp.getCurrentSlice() + "]");
    results.begin();
    // Perform the Gaussian fit
    long ellapsed = 0;
    if (!singleFit) {
        if (isLogProgress())
            IJ.log("Combined fit");
        // Estimate height from smoothed data
        double[] estimatedHeights = new double[maxIndices.length];
        for (int i = 0; i < estimatedHeights.length; i++) estimatedHeights[i] = smoothData[maxIndices[i]];
        FitConfiguration config = new FitConfiguration();
        setupPeakFiltering(config);
        long time = System.nanoTime();
        double[] params = fitMultiple(data, width, height, maxIndices, estimatedHeights);
        ellapsed = System.nanoTime() - time;
        if (params != null) {
            // Copy all the valid parameters into a new array
            double[] validParams = new double[params.length];
            int c = 0;
            int validPeaks = 0;
            validParams[c++] = params[0];
            double[] initialParams = convertParameters(fitResult.getInitialParameters());
            double[] paramsDev = convertParameters(fitResult.getParameterStdDev());
            Rectangle regionBounds = new Rectangle();
            int[] xpoints = new int[maxIndices.length];
            int[] ypoints = new int[maxIndices.length];
            int nMaxima = 0;
            for (int i = 1, n = 0; i < params.length; i += 6, n++) {
                int y = maxIndices[n] / width;
                int x = maxIndices[n] % width;
                // Check the peak is a good fit
                if (filterResults && config.validatePeak(n, initialParams, params) != FitStatus.OK)
                    continue;
                if (showFit) {
                    // Copy the valid parameters
                    validPeaks++;
                    for (int ii = i, j = 0; j < 6; ii++, j++) validParams[c++] = params[ii];
                }
                double[] peakParams = extractParams(params, i);
                double[] peakParamsDev = extractParams(paramsDev, i);
                addResult(bounds, regionBounds, data, peakParams, peakParamsDev, nMaxima, x, y, data[maxIndices[n]]);
                // Add fit result to the overlay - Coords are updated with the region offsets in addResult
                double xf = peakParams[3];
                double yf = peakParams[4];
                xpoints[nMaxima] = (int) (xf + 0.5);
                ypoints[nMaxima] = (int) (yf + 0.5);
                nMaxima++;
            }
            setOverlay(nMaxima, xpoints, ypoints);
            // Draw the fit
            if (showFit && validPeaks != 0) {
                double[] pixels = new double[data.length];
                EllipticalGaussian2DFunction f = new EllipticalGaussian2DFunction(validPeaks, width, height);
                invertParameters(validParams);
                f.initialise(validParams);
                for (int x = 0; x < pixels.length; x++) pixels[x] = f.eval(x);
                FloatProcessor fp = new FloatProcessor(width, height, pixels);
                // Insert into a full size image
                FloatProcessor fp2 = new FloatProcessor(ip.getWidth(), ip.getHeight());
                fp2.insert(fp, bounds.x, bounds.y);
                Utils.display(TITLE, fp2);
            }
        } else {
            if (isLogProgress()) {
                IJ.log("Failed to fit " + Utils.pleural(maxIndices.length, "peak") + getReason(fitResult));
            }
            imp.setOverlay(null);
        }
    } else {
        if (isLogProgress())
            IJ.log("Individual fit");
        int nMaxima = 0;
        int[] xpoints = new int[maxIndices.length];
        int[] ypoints = new int[maxIndices.length];
        // Extract each peak and fit individually
        ImageExtractor ie = new ImageExtractor(data, width, height);
        float[] region = null;
        Gaussian2DFitter gf = createGaussianFitter(filterResults);
        for (int n = 0; n < maxIndices.length; n++) {
            int y = maxIndices[n] / width;
            int x = maxIndices[n] % width;
            long time = System.nanoTime();
            Rectangle regionBounds = ie.getBoxRegionBounds(x, y, singleRegionSize);
            region = ie.crop(regionBounds, region);
            int newIndex = (y - regionBounds.y) * regionBounds.width + x - regionBounds.x;
            if (isLogProgress()) {
                IJ.log("Fitting peak " + (n + 1));
            }
            double[] peakParams = fitSingle(gf, region, regionBounds.width, regionBounds.height, newIndex, smoothData[maxIndices[n]]);
            ellapsed += System.nanoTime() - time;
            // Output fit result
            if (peakParams != null) {
                double[] peakParamsDev = null;
                if (showDeviations) {
                    peakParamsDev = convertParameters(fitResult.getParameterStdDev());
                }
                addResult(bounds, regionBounds, data, peakParams, peakParamsDev, n, x, y, data[maxIndices[n]]);
                // Add fit result to the overlay - Coords are updated with the region offsets in addResult
                double xf = peakParams[3];
                double yf = peakParams[4];
                xpoints[nMaxima] = (int) (xf + 0.5);
                ypoints[nMaxima] = (int) (yf + 0.5);
                nMaxima++;
            } else {
                if (isLogProgress()) {
                    IJ.log("Failed to fit peak " + (n + 1) + getReason(fitResult));
                }
            }
        }
        // Update the overlay
        if (nMaxima > 0)
            setOverlay(nMaxima, xpoints, ypoints);
        else
            imp.setOverlay(null);
    }
    results.end();
    if (isLogProgress())
        IJ.log("Time = " + (ellapsed / 1000000.0) + "ms");
}
Also used : FloatProcessor(ij.process.FloatProcessor) EllipticalGaussian2DFunction(gdsc.smlm.function.gaussian.EllipticalGaussian2DFunction) Gaussian2DFitter(gdsc.smlm.fitting.Gaussian2DFitter) Rectangle(java.awt.Rectangle) AverageFilter(gdsc.smlm.filters.AverageFilter) IJTablePeakResults(gdsc.smlm.ij.results.IJTablePeakResults) FitConfiguration(gdsc.smlm.fitting.FitConfiguration) GenericDialog(ij.gui.GenericDialog) ImageExtractor(gdsc.core.utils.ImageExtractor)

Example 2 with IJTablePeakResults

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

the class ResultsManager method addTableResults.

private void addTableResults(MemoryPeakResults results, PeakResultsList resultsList, boolean showDeviations, boolean showEndFrame) {
    if (resultsSettings.getResultsTable() != ResultsTable.NONE) {
        IJTablePeakResults r = new IJTablePeakResults(showDeviations);
        r.setPeakIdColumnName("Frame");
        r.setShowCalibratedValues(resultsSettings.getResultsTable() == ResultsTable.CALIBRATED);
        // Get a bias if required
        Calibration calibration = results.getCalibration();
        if (r.isShowCalibratedValues() && calibration.getBias() == 0) {
            GenericDialog gd = new GenericDialog(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()));
            }
        }
        r.setShowEndFrame(showEndFrame);
        resultsList.addOutput(r);
    }
}
Also used : IJTablePeakResults(gdsc.smlm.ij.results.IJTablePeakResults) ExtendedGenericDialog(ij.gui.ExtendedGenericDialog) GenericDialog(ij.gui.GenericDialog) Calibration(gdsc.smlm.results.Calibration)

Example 3 with IJTablePeakResults

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

the class PeakFit method addTableResults.

private void addTableResults(PeakResultsList resultsList) {
    if (resultsSettings.getResultsTable() != ResultsTable.NONE) {
        // imp.getTitle()
        String title = null;
        IJTablePeakResults r = new IJTablePeakResults(resultsSettings.showDeviations, title);
        r.setPeakIdColumnName("Frame");
        r.setShowCalibratedValues(resultsSettings.getResultsTable() == ResultsTable.CALIBRATED);
        r.setClearAtStart(simpleFit);
        r.setShowEndFrame(integrateFrames > 1);
        resultsList.addOutput(r);
    }
}
Also used : IJTablePeakResults(gdsc.smlm.ij.results.IJTablePeakResults)

Example 4 with IJTablePeakResults

use of gdsc.smlm.ij.results.IJTablePeakResults 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);
}
Also used : DescriptiveStatistics(org.apache.commons.math3.stat.descriptive.DescriptiveStatistics) ImageStack(ij.ImageStack) FloatProcessor(ij.process.FloatProcessor) Rectangle(java.awt.Rectangle) Calibration(gdsc.smlm.results.Calibration) ExtendedGenericDialog(ij.gui.ExtendedGenericDialog) ImagePlus(ij.ImagePlus) PeakResult(gdsc.smlm.results.PeakResult) ImageProcessor(ij.process.ImageProcessor) IJTablePeakResults(gdsc.smlm.ij.results.IJTablePeakResults) ImageSource(gdsc.smlm.results.ImageSource)

Aggregations

IJTablePeakResults (gdsc.smlm.ij.results.IJTablePeakResults)4 Calibration (gdsc.smlm.results.Calibration)2 ExtendedGenericDialog (ij.gui.ExtendedGenericDialog)2 GenericDialog (ij.gui.GenericDialog)2 FloatProcessor (ij.process.FloatProcessor)2 Rectangle (java.awt.Rectangle)2 ImageExtractor (gdsc.core.utils.ImageExtractor)1 AverageFilter (gdsc.smlm.filters.AverageFilter)1 FitConfiguration (gdsc.smlm.fitting.FitConfiguration)1 Gaussian2DFitter (gdsc.smlm.fitting.Gaussian2DFitter)1 EllipticalGaussian2DFunction (gdsc.smlm.function.gaussian.EllipticalGaussian2DFunction)1 ImageSource (gdsc.smlm.results.ImageSource)1 PeakResult (gdsc.smlm.results.PeakResult)1 ImagePlus (ij.ImagePlus)1 ImageStack (ij.ImageStack)1 ImageProcessor (ij.process.ImageProcessor)1 DescriptiveStatistics (org.apache.commons.math3.stat.descriptive.DescriptiveStatistics)1