Search in sources :

Example 31 with Plot

use of ij.gui.Plot in project GDSC-SMLM by aherbert.

the class DarkTimeAnalysis method analyse.

private void analyse(MemoryPeakResults results) {
    // Find min and max time frames
    results.sort();
    final int min = results.getFirstFrame();
    final int max = results.getLastFrame();
    // Trace results:
    // TODO - The search distance could have units to avoid assuming the results are in pixels
    final double d = settings.searchDistance / results.getCalibrationReader().getNmPerPixel();
    int range = max - min + 1;
    if (settings.maxDarkTime > 0) {
        range = Math.max(1, (int) Math.round(settings.maxDarkTime * 1000 / msPerFrame));
    }
    final TrackProgress tracker = SimpleImageJTrackProgress.getInstance();
    tracker.status("Analysing ...");
    tracker.log("Analysing (d=%s nm (%s px) t=%s s (%d frames)) ...", MathUtils.rounded(settings.searchDistance), MathUtils.rounded(d), MathUtils.rounded(range * msPerFrame / 1000.0), range);
    Trace[] traces;
    if (settings.method == 0) {
        final TraceManager tm = new TraceManager(results);
        tm.setTracker(tracker);
        tm.traceMolecules(d, range);
        traces = tm.getTraces();
    } else {
        final ClusteringEngine engine = new ClusteringEngine(Prefs.getThreads(), algorithms[settings.method - 1], tracker);
        final List<Cluster> clusters = engine.findClusters(TraceMolecules.convertToClusterPoints(results), d, range);
        traces = TraceMolecules.convertToTraces(results, clusters);
    }
    tracker.status("Computing histogram ...");
    // Build dark-time histogram
    final int[] times = new int[range];
    final StoredData stats = new StoredData();
    for (final Trace trace : traces) {
        if (trace.getBlinks() > 1) {
            for (final int t : trace.getOffTimes()) {
                times[t]++;
            }
            stats.add(trace.getOffTimes());
        }
    }
    plotDarkTimeHistogram(stats);
    // Cumulative histogram
    for (int i = 1; i < times.length; i++) {
        times[i] += times[i - 1];
    }
    final int total = times[times.length - 1];
    // Plot dark-time up to 100%
    double[] x = new double[range];
    double[] y = new double[range];
    int truncate = 0;
    for (int i = 0; i < x.length; i++) {
        x[i] = i * msPerFrame;
        if (times[i] == total) {
            // Final value at 100%
            y[i] = 100.0;
            truncate = i + 1;
            break;
        }
        y[i] = (100.0 * times[i]) / total;
    }
    if (truncate > 0) {
        x = Arrays.copyOf(x, truncate);
        y = Arrays.copyOf(y, truncate);
    }
    final String title = "Cumulative Dark-time";
    final Plot plot = new Plot(title, "Time (ms)", "Percentile");
    plot.addPoints(x, y, Plot.LINE);
    ImageJUtils.display(title, plot);
    // Report percentile
    for (int i = 0; i < y.length; i++) {
        if (y[i] >= settings.percentile) {
            ImageJUtils.log("Dark-time Percentile %.1f @ %s ms = %s s", settings.percentile, MathUtils.rounded(x[i]), MathUtils.rounded(x[i] / 1000));
            break;
        }
    }
    tracker.status("");
}
Also used : Plot(ij.gui.Plot) SimpleImageJTrackProgress(uk.ac.sussex.gdsc.core.ij.SimpleImageJTrackProgress) TrackProgress(uk.ac.sussex.gdsc.core.logging.TrackProgress) Cluster(uk.ac.sussex.gdsc.core.clustering.Cluster) TraceManager(uk.ac.sussex.gdsc.smlm.results.TraceManager) Trace(uk.ac.sussex.gdsc.smlm.results.Trace) StoredData(uk.ac.sussex.gdsc.core.utils.StoredData) ClusteringEngine(uk.ac.sussex.gdsc.core.clustering.ClusteringEngine)

Example 32 with Plot

use of ij.gui.Plot in project GDSC-SMLM by aherbert.

the class AstigmatismModelManager method plot.

private static PlotWindow plot(WindowOrganiser wo, double[] z, String yTitle, double[] y1, String y1Title, double[] y2, String y2Title) {
    final String title = TITLE + " " + yTitle;
    final Plot plot = new Plot(title, "Z (μm)", yTitle);
    double[] limits = MathUtils.limits(y1);
    if (y2 != null) {
        limits = MathUtils.limits(limits, y2);
    }
    final double rangex = (z[z.length - 1] - z[0]) * 0.05;
    final double rangey = (limits[1] - limits[0]) * 0.05;
    plot.setLimits(z[0] - rangex, z[z.length - 1] + rangex, limits[0] - rangey, limits[1] + rangey);
    plot.setColor(Color.RED);
    plot.addPoints(z, y1, Plot.CIRCLE);
    if (y2 != null) {
        plot.setColor(Color.BLUE);
        plot.addPoints(z, y2, Plot.CIRCLE);
        plot.setColor(Color.BLACK);
        plot.addLegend(y1Title + "\n" + y2Title);
    }
    return ImageJUtils.display(title, plot, 0, wo);
}
Also used : Plot(ij.gui.Plot)

Example 33 with Plot

use of ij.gui.Plot in project GDSC-SMLM by aherbert.

the class AstigmatismModelManager method viewModel.

private void viewModel() {
    final ExtendedGenericDialog gd = new ExtendedGenericDialog(TITLE);
    final String[] models = listAstigmatismModels(false);
    gd.addChoice("Model", models, pluginSettings.getSelected());
    gd.addChoice("z_distance_unit", SettingsManager.getDistanceUnitNames(), pluginSettings.getZDistanceUnitValue());
    gd.addChoice("s_distance_unit", SettingsManager.getDistanceUnitNames(), pluginSettings.getSDistanceUnitValue());
    gd.addCheckbox("Show_depth_of_focus", pluginSettings.getShowDepthOfFocus());
    gd.addCheckbox("Show_combined_width", pluginSettings.getShowCombinedWidth());
    gd.addCheckbox("Show_PSF", pluginSettings.getShowPsf());
    gd.addHelp(HelpUrls.getUrl("astigmatism-model-manager-view"));
    gd.showDialog();
    if (gd.wasCanceled()) {
        return;
    }
    final String name = gd.getNextChoice();
    pluginSettings.setSelected(name);
    pluginSettings.setZDistanceUnitValue(gd.getNextChoiceIndex());
    pluginSettings.setSDistanceUnitValue(gd.getNextChoiceIndex());
    pluginSettings.setShowDepthOfFocus(gd.getNextBoolean());
    pluginSettings.setShowCombinedWidth(gd.getNextBoolean());
    pluginSettings.setShowPsf(gd.getNextBoolean());
    // Try and get the named resource
    AstigmatismModel model = AstigmatismModelSettingsHolder.getSettings().getAstigmatismModelResourcesMap().get(name);
    if (model == null) {
        IJ.error(TITLE, "Failed to find astigmatism model: " + name);
        return;
    }
    try {
        model = convert(model, pluginSettings.getZDistanceUnit(), pluginSettings.getSDistanceUnit());
    } catch (final ConversionException ex) {
        ImageJUtils.log("Bad conversion (%s), defaulting to native model units", ex.getMessage());
    }
    ImageJUtils.log("Astigmatism model: %s\n%s", name, model);
    // Plot the curve. Do this so we encompass twice the depth-of-field.
    final double gamma = model.getGamma();
    final double d = model.getD();
    final double s0x = model.getS0X();
    final double Ax = model.getAx();
    final double Bx = model.getBx();
    final double s0y = model.getS0Y();
    final double Ay = model.getAy();
    final double By = model.getBy();
    final double range = Math.abs(gamma) + 1.5 * d;
    final int n = 200;
    final double step = range / n;
    final double[] z = new double[2 * n + 1];
    final double[] sx = new double[z.length];
    final double[] sy = new double[z.length];
    // Use the same class that is used during fitting
    final HoltzerAstigmatismZModel m = HoltzerAstigmatismZModel.create(s0x, s0y, gamma, d, Ax, Bx, Ay, By);
    for (int i = 0; i < z.length; i++) {
        final double zz = -range + i * step;
        z[i] = zz;
        sx[i] = m.getSx(zz);
        sy[i] = m.getSy(zz);
    }
    final String title = TITLE + " Width Curve";
    final Plot plot = new Plot(title, "Z (" + UnitHelper.getShortName(model.getZDistanceUnit()) + ")", "Width (" + UnitHelper.getShortName(model.getSDistanceUnit()) + ")");
    double[] limits = MathUtils.limits(sx);
    limits = MathUtils.limits(limits, sy);
    final double rangex = (z[z.length - 1] - z[0]) * 0.05;
    final double rangey = (limits[1] - limits[0]) * 0.05;
    final double miny = limits[0] - rangey;
    final double maxy = limits[1] + rangey;
    plot.setLimits(z[0] - rangex, z[z.length - 1] + rangex, miny, maxy);
    plot.setColor(Color.RED);
    plot.addPoints(z, sx, Plot.LINE);
    plot.setColor(Color.BLUE);
    plot.addPoints(z, sy, Plot.LINE);
    plot.setColor(Color.YELLOW);
    if (pluginSettings.getShowDepthOfFocus()) {
        final double z0x = gamma;
        final double z0y = -gamma;
        plot.setColor(Color.RED.darker());
        plot.drawDottedLine(z0x - d, miny, z0x - d, maxy, 4);
        plot.drawDottedLine(z0x + d, miny, z0x + d, maxy, 4);
        plot.setColor(Color.BLUE.darker());
        plot.drawDottedLine(z0y - d, miny, z0y - d, maxy, 4);
        plot.drawDottedLine(z0y + d, miny, z0y + d, maxy, 4);
    }
    String legend = "Sx\nSy";
    if (pluginSettings.getShowCombinedWidth()) {
        final double[] s = new double[z.length];
        for (int i = 0; i < z.length; i++) {
            s[i] = Gaussian2DPeakResultHelper.getStandardDeviation(sx[i], sy[i]);
        }
        plot.setColor(Color.GREEN);
        plot.addPoints(z, s, Plot.LINE);
        legend += "\tS";
    }
    plot.setColor(Color.BLACK);
    plot.addLegend(legend);
    plot.addLabel(0, 0, String.format("Model = %s (%s nm/pixel)", name, MathUtils.rounded(model.getNmPerPixel())));
    ImageJUtils.display(title, plot);
    if (!pluginSettings.getShowPsf()) {
        return;
    }
    // Get pixel range using 3x[max SD]
    final int width = 1 + 2 * ((int) Math.ceil(limits[1] * 3));
    new ModelRenderer(name, model, m, range, width, plot).run();
}
Also used : ConversionException(uk.ac.sussex.gdsc.core.data.utils.ConversionException) AstigmatismModel(uk.ac.sussex.gdsc.smlm.data.config.PSFProtos.AstigmatismModel) Plot(ij.gui.Plot) NonBlockingExtendedGenericDialog(uk.ac.sussex.gdsc.core.ij.gui.NonBlockingExtendedGenericDialog) ExtendedGenericDialog(uk.ac.sussex.gdsc.core.ij.gui.ExtendedGenericDialog) HoltzerAstigmatismZModel(uk.ac.sussex.gdsc.smlm.function.gaussian.HoltzerAstigmatismZModel)

Example 34 with Plot

use of ij.gui.Plot in project GDSC-SMLM by aherbert.

the class AstigmatismModelManager method doCurveFit.

private boolean doCurveFit() {
    // Estimate:
    // Focal plane = where width is at a minimum
    // s0x/s0y = the min width of x/y
    // gamma = Half the distance between the focal planes
    // z0 = half way between the two focal planes
    // d = depth of focus
    double[] smoothSx = fitSx;
    double[] smoothSy = fitSy;
    if (pluginSettings.getSmoothing() > 0) {
        final LoessInterpolator loess = new LoessInterpolator(pluginSettings.getSmoothing(), 0);
        smoothSx = loess.smooth(fitZ, fitSx);
        smoothSy = loess.smooth(fitZ, fitSy);
        final Plot plot = widthPlot.getPlot();
        plot.setColor(Color.RED);
        plot.addPoints(fitZ, smoothSx, Plot.LINE);
        plot.setColor(Color.BLUE);
        plot.addPoints(fitZ, smoothSy, Plot.LINE);
        plot.setColor(Color.BLACK);
        plot.updateImage();
    }
    final int focalPlaneXindex = SimpleArrayUtils.findMinIndex(smoothSx);
    final int focalPlaneYindex = SimpleArrayUtils.findMinIndex(smoothSy);
    final double s0x = smoothSx[focalPlaneXindex];
    final double s0y = smoothSy[focalPlaneYindex];
    final double focalPlaneX = fitZ[focalPlaneXindex];
    final double focalPlaneY = fitZ[focalPlaneYindex];
    double gamma = Math.abs(focalPlaneY - focalPlaneX) / 2;
    final double z0 = (focalPlaneX + focalPlaneY) / 2;
    final double d = (estimateD(focalPlaneXindex, fitZ, smoothSx) + estimateD(focalPlaneYindex, fitZ, smoothSy)) / 2;
    // Start with Ax, Bx, Ay, By as zero.
    final double Ax = 0;
    final double Bx = 0;
    final double Ay = 0;
    final double By = 0;
    // If this is not the case we can invert the gamma parameter.
    if (focalPlaneXindex < focalPlaneYindex) {
        gamma = -gamma;
    }
    // Use an LVM fitter with numerical gradients.
    final double initialStepBoundFactor = 100;
    final double costRelativeTolerance = 1e-10;
    final double parRelativeTolerance = 1e-10;
    final double orthoTolerance = 1e-10;
    final double threshold = Precision.SAFE_MIN;
    // We optimise against both sx and sy as a combined y-value.
    final double[] y = new double[fitZ.length * 2];
    System.arraycopy(fitSx, 0, y, 0, fitSx.length);
    System.arraycopy(fitSy, 0, y, fitSx.length, fitSy.length);
    final LevenbergMarquardtOptimizer optimizer = new LevenbergMarquardtOptimizer(initialStepBoundFactor, costRelativeTolerance, parRelativeTolerance, orthoTolerance, threshold);
    parameters = new double[9];
    parameters[P_GAMMA] = gamma;
    parameters[P_Z0] = z0;
    parameters[P_D] = d;
    parameters[P_S0X] = s0x;
    parameters[P_AX] = Ax;
    parameters[P_BX] = Bx;
    parameters[P_S0Y] = s0y;
    parameters[P_AY] = Ay;
    parameters[P_BY] = By;
    record("Initial", parameters);
    if (pluginSettings.getShowEstimatedCurve()) {
        plotFit(parameters);
        IJ.showMessage(TITLE, "Showing the estimated curve parameters.\nClick OK to continue.");
    }
    // @formatter:off
    final LeastSquaresBuilder builder = new LeastSquaresBuilder().maxEvaluations(Integer.MAX_VALUE).maxIterations(3000).start(parameters).target(y);
    if (pluginSettings.getWeightedFit()) {
        builder.weight(new DiagonalMatrix(getWeights(smoothSx, smoothSy)));
    }
    final AstigmatismVectorFunction vf = new AstigmatismVectorFunction();
    builder.model(vf, new AstigmatismMatrixFunction());
    final LeastSquaresProblem problem = builder.build();
    try {
        final Optimum optimum = optimizer.optimize(problem);
        parameters = optimum.getPoint().toArray();
        record("Final", parameters);
        plotFit(parameters);
        saveResult(optimum);
    } catch (final Exception ex) {
        IJ.error(TITLE, "Failed to fit curve: " + ex.getMessage());
        return false;
    }
    return true;
}
Also used : Plot(ij.gui.Plot) ConversionException(uk.ac.sussex.gdsc.core.data.utils.ConversionException) LeastSquaresBuilder(org.apache.commons.math3.fitting.leastsquares.LeastSquaresBuilder) LoessInterpolator(org.apache.commons.math3.analysis.interpolation.LoessInterpolator) Optimum(org.apache.commons.math3.fitting.leastsquares.LeastSquaresOptimizer.Optimum) LevenbergMarquardtOptimizer(org.apache.commons.math3.fitting.leastsquares.LevenbergMarquardtOptimizer) DiagonalMatrix(org.apache.commons.math3.linear.DiagonalMatrix) LeastSquaresProblem(org.apache.commons.math3.fitting.leastsquares.LeastSquaresProblem)

Example 35 with Plot

use of ij.gui.Plot in project GDSC-SMLM by aherbert.

the class CameraModelFisherInformationAnalysis method plotFromCache.

private static void plotFromCache() {
    // Build a list of curve stored in the cache
    final String[] names = new String[cache.size()];
    final PoissonFisherInformationData[] datas = new PoissonFisherInformationData[cache.size()];
    int count = 0;
    for (final Entry<FiKey, PoissonFisherInformationData> e : cache.entrySet()) {
        names[count] = e.getKey().toString();
        datas[count] = e.getValue();
        count++;
    }
    final ExtendedGenericDialog gd = new ExtendedGenericDialog(TITLE);
    gd.addChoice("Fisher_information", names, cachePlot);
    gd.addChoice("Plot_point", POINT_OPTION, pointOption);
    gd.showDialog();
    if (gd.wasCanceled()) {
        return;
    }
    final int index = gd.getNextChoiceIndex();
    pointOption = gd.getNextChoiceIndex();
    cachePlot = names[index];
    final PoissonFisherInformationData data = datas[index];
    final FiKey key = new FiKey(data);
    count = 0;
    final double[] exp = new double[data.getAlphaSampleCount()];
    final double[] alpha = new double[data.getAlphaSampleCount()];
    for (final AlphaSample s : data.getAlphaSampleList()) {
        exp[count] = s.getLog10Mean();
        alpha[count] = s.getAlpha();
        count++;
    }
    // Just in case
    // Sort.sortArrays(alpha, exp, true);
    // Test if we can use ImageJ support for a X log scale
    final boolean logScaleX = ((float) Math.pow(10, exp[0]) != 0);
    double[] x = exp;
    String xTitle = "log10(photons)";
    if (logScaleX) {
        final double[] photons = new double[exp.length];
        for (int i = 0; i < photons.length; i++) {
            photons[i] = Math.pow(10, exp[i]);
        }
        x = photons;
        xTitle = "photons";
    }
    // Get interpolation for alpha. Convert to base e.
    final double[] logU = exp.clone();
    final double scale = Math.log(10);
    for (int i = 0; i < logU.length; i++) {
        logU[i] *= scale;
    }
    final BasePoissonFisherInformation if1 = getInterpolatedPoissonFisherInformation(key.getType(), logU, alpha, null);
    // Interpolate with 5 points per sample for smooth curve
    final int n = 5;
    final TDoubleArrayList iexp = new TDoubleArrayList();
    final TDoubleArrayList iphotons = new TDoubleArrayList();
    for (int i = 1; i < exp.length; i++) {
        final int i_1 = i - 1;
        final double h = (exp[i] - exp[i_1]) / n;
        for (int j = 0; j < n; j++) {
            final double e = exp[i_1] + j * h;
            iexp.add(e);
            iphotons.add(Math.pow(10, e));
        }
    }
    iexp.add(exp[exp.length - 1]);
    iphotons.add(Math.pow(10, exp[exp.length - 1]));
    final double[] photons = iphotons.toArray();
    final double[] ix = (logScaleX) ? photons : iexp.toArray();
    final double[] ialpha1 = getAlpha(if1, photons);
    final int pointShape = getPointShape(pointOption);
    final String title = "Cached Relative Fisher Information";
    final Plot plot = new Plot(title, xTitle, "Noise coefficient (alpha)");
    plot.setLimits(x[0], x[x.length - 1], -0.05, 1.05);
    if (logScaleX) {
        plot.setLogScaleX();
    }
    plot.setColor(Color.blue);
    plot.addPoints(ix, ialpha1, Plot.LINE);
    // Option to show nodes
    if (pointShape != -1) {
        plot.addPoints(x, alpha, pointShape);
    }
    plot.setColor(Color.BLACK);
    plot.addLabel(0, 0, cachePlot);
    ImageJUtils.display(title, plot);
}
Also used : Plot(ij.gui.Plot) AlphaSample(uk.ac.sussex.gdsc.smlm.data.config.FisherProtos.AlphaSample) NonBlockingExtendedGenericDialog(uk.ac.sussex.gdsc.core.ij.gui.NonBlockingExtendedGenericDialog) ExtendedGenericDialog(uk.ac.sussex.gdsc.core.ij.gui.ExtendedGenericDialog) PoissonFisherInformationData(uk.ac.sussex.gdsc.smlm.data.config.FisherProtos.PoissonFisherInformationData) TDoubleArrayList(gnu.trove.list.array.TDoubleArrayList) BasePoissonFisherInformation(uk.ac.sussex.gdsc.smlm.function.BasePoissonFisherInformation)

Aggregations

Plot (ij.gui.Plot)89 HistogramPlot (uk.ac.sussex.gdsc.core.ij.HistogramPlot)20 Point (java.awt.Point)19 PlotWindow (ij.gui.PlotWindow)17 Color (java.awt.Color)13 WindowOrganiser (uk.ac.sussex.gdsc.core.ij.plugin.WindowOrganiser)13 HistogramPlotBuilder (uk.ac.sussex.gdsc.core.ij.HistogramPlot.HistogramPlotBuilder)12 BasePoint (uk.ac.sussex.gdsc.core.match.BasePoint)12 ExtendedGenericDialog (uk.ac.sussex.gdsc.core.ij.gui.ExtendedGenericDialog)11 MemoryPeakResults (uk.ac.sussex.gdsc.smlm.results.MemoryPeakResults)11 Rectangle (java.awt.Rectangle)9 ArrayList (java.util.ArrayList)9 GenericDialog (ij.gui.GenericDialog)8 NonBlockingExtendedGenericDialog (uk.ac.sussex.gdsc.core.ij.gui.NonBlockingExtendedGenericDialog)7 LocalList (uk.ac.sussex.gdsc.core.utils.LocalList)7 Statistics (uk.ac.sussex.gdsc.core.utils.Statistics)7 StoredData (uk.ac.sussex.gdsc.core.utils.StoredData)7 StoredDataStatistics (uk.ac.sussex.gdsc.core.utils.StoredDataStatistics)7 ImagePlus (ij.ImagePlus)6 TDoubleArrayList (gnu.trove.list.array.TDoubleArrayList)5