Search in sources :

Example 1 with WindowOrganiser

use of ij.plugin.WindowOrganiser 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;
}
Also used : InterlacedImageSource(gdsc.smlm.results.InterlacedImageSource) AggregatedImageSource(gdsc.smlm.results.AggregatedImageSource) ImageStack(ij.ImageStack) Rectangle(java.awt.Rectangle) StoredDataStatistics(gdsc.core.utils.StoredDataStatistics) WindowOrganiser(ij.plugin.WindowOrganiser) IJImageSource(gdsc.smlm.ij.IJImageSource) ImageProcessor(ij.process.ImageProcessor) FitEngine(gdsc.smlm.engine.FitEngine) Random(gdsc.core.utils.Random) FitConfiguration(gdsc.smlm.fitting.FitConfiguration) InterlacedImageSource(gdsc.smlm.results.InterlacedImageSource) ImageSource(gdsc.smlm.results.ImageSource) AggregatedImageSource(gdsc.smlm.results.AggregatedImageSource) IJImageSource(gdsc.smlm.ij.IJImageSource) FitJob(gdsc.smlm.engine.FitJob)

Example 2 with WindowOrganiser

use of ij.plugin.WindowOrganiser in project GDSC-SMLM by aherbert.

the class FilterAnalysis method showPlots.

private void showPlots() {
    if (plots.isEmpty())
        return;
    // Display the top N plots
    int[] list = new int[plots.size()];
    int i = 0;
    for (NamedPlot p : plots) {
        Plot2 plot = new Plot2(p.name, p.xAxisName, "Jaccard", p.xValues, p.yValues);
        plot.setLimits(p.xValues[0], p.xValues[p.xValues.length - 1], 0, 1);
        plot.setColor(Color.RED);
        plot.draw();
        plot.setColor(Color.BLUE);
        plot.addPoints(p.xValues, p.yValues, Plot2.CROSS);
        PlotWindow plotWindow = Utils.display(p.name, plot);
        list[i++] = plotWindow.getImagePlus().getID();
    }
    new WindowOrganiser().tileWindows(list);
}
Also used : PlotWindow(ij.gui.PlotWindow) Plot2(ij.gui.Plot2) WindowOrganiser(ij.plugin.WindowOrganiser)

Example 3 with WindowOrganiser

use of ij.plugin.WindowOrganiser in project GDSC-SMLM by aherbert.

the class DiffusionRateTest method run.

/*
	 * (non-Javadoc)
	 * 
	 * @see ij.plugin.PlugIn#run(java.lang.String)
	 */
public void run(String arg) {
    SMLMUsageTracker.recordPlugin(this.getClass(), arg);
    if (IJ.controlKeyDown()) {
        simpleTest();
        return;
    }
    extraOptions = Utils.isExtraOptions();
    if (!showDialog())
        return;
    lastSimulatedDataset[0] = lastSimulatedDataset[1] = "";
    lastSimulatedPrecision = 0;
    final int totalSteps = (int) Math.ceil(settings.seconds * settings.stepsPerSecond);
    conversionFactor = 1000000.0 / (settings.pixelPitch * settings.pixelPitch);
    // Diffusion rate is um^2/sec. Convert to pixels per simulation frame.
    final double diffusionRateInPixelsPerSecond = settings.diffusionRate * conversionFactor;
    final double diffusionRateInPixelsPerStep = diffusionRateInPixelsPerSecond / settings.stepsPerSecond;
    final double precisionInPixels = myPrecision / settings.pixelPitch;
    final boolean addError = myPrecision != 0;
    Utils.log(TITLE + " : D = %s um^2/sec, Precision = %s nm", Utils.rounded(settings.diffusionRate, 4), Utils.rounded(myPrecision, 4));
    Utils.log("Mean-displacement per dimension = %s nm/sec", Utils.rounded(1e3 * ImageModel.getRandomMoveDistance(settings.diffusionRate), 4));
    if (extraOptions)
        Utils.log("Step size = %s, precision = %s", Utils.rounded(ImageModel.getRandomMoveDistance(diffusionRateInPixelsPerStep)), Utils.rounded(precisionInPixels));
    // Convert diffusion co-efficient into the standard deviation for the random walk
    final double diffusionSigma = (settings.getDiffusionType() == DiffusionType.LINEAR_WALK) ? // Q. What should this be? At the moment just do 1D diffusion on a random vector
    ImageModel.getRandomMoveDistance(diffusionRateInPixelsPerStep) : ImageModel.getRandomMoveDistance(diffusionRateInPixelsPerStep);
    Utils.log("Simulation step-size = %s nm", Utils.rounded(settings.pixelPitch * diffusionSigma, 4));
    // Move the molecules and get the diffusion rate
    IJ.showStatus("Simulating ...");
    final long start = System.nanoTime();
    final long seed = System.currentTimeMillis() + System.identityHashCode(this);
    RandomGenerator[] random = new RandomGenerator[3];
    RandomGenerator[] random2 = new RandomGenerator[3];
    for (int i = 0; i < 3; i++) {
        random[i] = new Well19937c(seed + i * 12436);
        random2[i] = new Well19937c(seed + i * 678678 + 3);
    }
    Statistics[] stats2D = new Statistics[totalSteps];
    Statistics[] stats3D = new Statistics[totalSteps];
    StoredDataStatistics jumpDistances2D = new StoredDataStatistics(totalSteps);
    StoredDataStatistics jumpDistances3D = new StoredDataStatistics(totalSteps);
    for (int j = 0; j < totalSteps; j++) {
        stats2D[j] = new Statistics();
        stats3D[j] = new Statistics();
    }
    SphericalDistribution dist = new SphericalDistribution(settings.confinementRadius / settings.pixelPitch);
    Statistics asymptote = new Statistics();
    // Save results to memory
    MemoryPeakResults results = new MemoryPeakResults(totalSteps);
    Calibration cal = new Calibration(settings.pixelPitch, 1, 1000.0 / settings.stepsPerSecond);
    results.setCalibration(cal);
    results.setName(TITLE);
    int peak = 0;
    // Store raw coordinates
    ArrayList<Point> points = new ArrayList<Point>(totalSteps);
    StoredData totalJumpDistances1D = new StoredData(settings.particles);
    StoredData totalJumpDistances2D = new StoredData(settings.particles);
    StoredData totalJumpDistances3D = new StoredData(settings.particles);
    for (int i = 0; i < settings.particles; i++) {
        if (i % 16 == 0) {
            IJ.showProgress(i, settings.particles);
            if (Utils.isInterrupted())
                return;
        }
        // Increment the frame so that tracing analysis can distinguish traces
        peak++;
        double[] origin = new double[3];
        final int id = i + 1;
        MoleculeModel m = new MoleculeModel(id, origin.clone());
        if (addError)
            origin = addError(origin, precisionInPixels, random);
        if (useConfinement) {
            // Note: When using confinement the average displacement should asymptote
            // at the average distance of a point from the centre of a ball. This is 3r/4.
            // See: http://answers.yahoo.com/question/index?qid=20090131162630AAMTUfM
            // The equivalent in 2D is 2r/3. However although we are plotting 2D distance
            // this is a projection of the 3D position onto the plane and so the particles
            // will not be evenly spread (there will be clustering at centre caused by the
            // poles)
            final double[] axis = (settings.getDiffusionType() == DiffusionType.LINEAR_WALK) ? nextVector() : null;
            for (int j = 0; j < totalSteps; j++) {
                double[] xyz = m.getCoordinates();
                double[] originalXyz = xyz.clone();
                for (int n = confinementAttempts; n-- > 0; ) {
                    if (settings.getDiffusionType() == DiffusionType.GRID_WALK)
                        m.walk(diffusionSigma, random);
                    else if (settings.getDiffusionType() == DiffusionType.LINEAR_WALK)
                        m.slide(diffusionSigma, axis, random[0]);
                    else
                        m.move(diffusionSigma, random);
                    if (!dist.isWithin(m.getCoordinates())) {
                        // Reset position
                        for (int k = 0; k < 3; k++) xyz[k] = originalXyz[k];
                    } else {
                        // The move was allowed
                        break;
                    }
                }
                points.add(new Point(id, xyz));
                if (addError)
                    xyz = addError(xyz, precisionInPixels, random2);
                peak = record(xyz, id, peak, stats2D[j], stats3D[j], jumpDistances2D, jumpDistances3D, origin, results);
            }
            asymptote.add(distance(m.getCoordinates()));
        } else {
            if (settings.getDiffusionType() == DiffusionType.GRID_WALK) {
                for (int j = 0; j < totalSteps; j++) {
                    m.walk(diffusionSigma, random);
                    double[] xyz = m.getCoordinates();
                    points.add(new Point(id, xyz));
                    if (addError)
                        xyz = addError(xyz, precisionInPixels, random2);
                    peak = record(xyz, id, peak, stats2D[j], stats3D[j], jumpDistances2D, jumpDistances3D, origin, results);
                }
            } else if (settings.getDiffusionType() == DiffusionType.LINEAR_WALK) {
                final double[] axis = nextVector();
                for (int j = 0; j < totalSteps; j++) {
                    m.slide(diffusionSigma, axis, random[0]);
                    double[] xyz = m.getCoordinates();
                    points.add(new Point(id, xyz));
                    if (addError)
                        xyz = addError(xyz, precisionInPixels, random2);
                    peak = record(xyz, id, peak, stats2D[j], stats3D[j], jumpDistances2D, jumpDistances3D, origin, results);
                }
            } else {
                for (int j = 0; j < totalSteps; j++) {
                    m.move(diffusionSigma, random);
                    double[] xyz = m.getCoordinates();
                    points.add(new Point(id, xyz));
                    if (addError)
                        xyz = addError(xyz, precisionInPixels, random2);
                    peak = record(xyz, id, peak, stats2D[j], stats3D[j], jumpDistances2D, jumpDistances3D, origin, results);
                }
            }
        }
        // Debug: record all the particles so they can be analysed
        // System.out.printf("%f %f %f\n", m.getX(), m.getY(), m.getZ());
        final double[] xyz = m.getCoordinates();
        double d2 = 0;
        totalJumpDistances1D.add(d2 = xyz[0] * xyz[0]);
        totalJumpDistances2D.add(d2 += xyz[1] * xyz[1]);
        totalJumpDistances3D.add(d2 += xyz[2] * xyz[2]);
    }
    final double time = (System.nanoTime() - start) / 1000000.0;
    IJ.showProgress(1);
    MemoryPeakResults.addResults(results);
    lastSimulatedDataset[0] = results.getName();
    lastSimulatedPrecision = myPrecision;
    // Convert pixels^2/step to um^2/sec
    final double msd2D = (jumpDistances2D.getMean() / conversionFactor) / (results.getCalibration().getExposureTime() / 1000);
    final double msd3D = (jumpDistances3D.getMean() / conversionFactor) / (results.getCalibration().getExposureTime() / 1000);
    Utils.log("Raw data D=%s um^2/s, Precision = %s nm, N=%d, step=%s s, mean2D=%s um^2, MSD 2D = %s um^2/s, mean3D=%s um^2, MSD 3D = %s um^2/s", Utils.rounded(settings.diffusionRate), Utils.rounded(myPrecision), jumpDistances2D.getN(), Utils.rounded(results.getCalibration().getExposureTime() / 1000), Utils.rounded(jumpDistances2D.getMean() / conversionFactor), Utils.rounded(msd2D), Utils.rounded(jumpDistances3D.getMean() / conversionFactor), Utils.rounded(msd3D));
    aggregateIntoFrames(points, addError, precisionInPixels, random2);
    IJ.showStatus("Analysing results ...");
    if (showDiffusionExample) {
        showExample(totalSteps, diffusionSigma, random);
    }
    // Plot a graph of mean squared distance
    double[] xValues = new double[stats2D.length];
    double[] yValues2D = new double[stats2D.length];
    double[] yValues3D = new double[stats3D.length];
    double[] upper2D = new double[stats2D.length];
    double[] lower2D = new double[stats2D.length];
    double[] upper3D = new double[stats3D.length];
    double[] lower3D = new double[stats3D.length];
    SimpleRegression r2D = new SimpleRegression(false);
    SimpleRegression r3D = new SimpleRegression(false);
    final int firstN = (useConfinement) ? fitN : totalSteps;
    for (int j = 0; j < totalSteps; j++) {
        // Convert steps to seconds
        xValues[j] = (double) (j + 1) / settings.stepsPerSecond;
        // Convert values in pixels^2 to um^2
        final double mean2D = stats2D[j].getMean() / conversionFactor;
        final double mean3D = stats3D[j].getMean() / conversionFactor;
        final double sd2D = stats2D[j].getStandardDeviation() / conversionFactor;
        final double sd3D = stats3D[j].getStandardDeviation() / conversionFactor;
        yValues2D[j] = mean2D;
        yValues3D[j] = mean3D;
        upper2D[j] = mean2D + sd2D;
        lower2D[j] = mean2D - sd2D;
        upper3D[j] = mean3D + sd3D;
        lower3D[j] = mean3D - sd3D;
        if (j < firstN) {
            r2D.addData(xValues[j], yValues2D[j]);
            r3D.addData(xValues[j], yValues3D[j]);
        }
    }
    // TODO - Fit using the equation for 2D confined diffusion:
    // MSD = 4s^2 + R^2 (1 - 0.99e^(-1.84^2 Dt / R^2)
    // s = localisation precision
    // R = confinement radius
    // D = 2D diffusion coefficient
    // t = time
    final PolynomialFunction fitted2D, fitted3D;
    if (r2D.getN() > 0) {
        // Do linear regression to get diffusion rate
        final double[] best2D = new double[] { r2D.getIntercept(), r2D.getSlope() };
        fitted2D = new PolynomialFunction(best2D);
        final double[] best3D = new double[] { r3D.getIntercept(), r3D.getSlope() };
        fitted3D = new PolynomialFunction(best3D);
        // For 2D diffusion: d^2 = 4D
        // where: d^2 = mean-square displacement
        double D = best2D[1] / 4.0;
        String msg = "2D Diffusion rate = " + Utils.rounded(D, 4) + " um^2 / sec (" + Utils.timeToString(time) + ")";
        IJ.showStatus(msg);
        Utils.log(msg);
        D = best3D[1] / 6.0;
        Utils.log("3D Diffusion rate = " + Utils.rounded(D, 4) + " um^2 / sec (" + Utils.timeToString(time) + ")");
    } else {
        fitted2D = fitted3D = null;
    }
    // Create plots
    plotMSD(totalSteps, xValues, yValues2D, lower2D, upper2D, fitted2D, 2);
    plotMSD(totalSteps, xValues, yValues3D, lower3D, upper3D, fitted3D, 3);
    plotJumpDistances(TITLE, jumpDistances2D, 2, 1);
    plotJumpDistances(TITLE, jumpDistances3D, 3, 1);
    if (idCount > 0)
        new WindowOrganiser().tileWindows(idList);
    if (useConfinement)
        Utils.log("3D asymptote distance = %s nm (expected %.2f)", Utils.rounded(asymptote.getMean() * settings.pixelPitch, 4), 3 * settings.confinementRadius / 4);
}
Also used : SphericalDistribution(gdsc.smlm.model.SphericalDistribution) StoredDataStatistics(gdsc.core.utils.StoredDataStatistics) ArrayList(java.util.ArrayList) PolynomialFunction(org.apache.commons.math3.analysis.polynomials.PolynomialFunction) Calibration(gdsc.smlm.results.Calibration) WindowOrganiser(ij.plugin.WindowOrganiser) Well19937c(org.apache.commons.math3.random.Well19937c) StoredDataStatistics(gdsc.core.utils.StoredDataStatistics) Statistics(gdsc.core.utils.Statistics) RandomGenerator(org.apache.commons.math3.random.RandomGenerator) MoleculeModel(gdsc.smlm.model.MoleculeModel) SimpleRegression(org.apache.commons.math3.stat.regression.SimpleRegression) StoredData(gdsc.core.utils.StoredData) MemoryPeakResults(gdsc.smlm.results.MemoryPeakResults)

Example 4 with WindowOrganiser

use of ij.plugin.WindowOrganiser in project GDSC-SMLM by aherbert.

the class DiffusionRateTest method simpleTest.

/**
	 * Perform a simple diffusion test. This can be used to understand the distributions that are generated during
	 * 3D diffusion.
	 */
private void simpleTest() {
    if (!showSimpleDialog())
        return;
    StoredDataStatistics[] stats2 = new StoredDataStatistics[3];
    StoredDataStatistics[] stats = new StoredDataStatistics[3];
    RandomGenerator[] random = new RandomGenerator[3];
    final long seed = System.currentTimeMillis() + System.identityHashCode(this);
    for (int i = 0; i < 3; i++) {
        stats2[i] = new StoredDataStatistics(simpleParticles);
        stats[i] = new StoredDataStatistics(simpleParticles);
        random[i] = new Well19937c(seed + i);
    }
    final double scale = Math.sqrt(2 * simpleD);
    final int report = Math.max(1, simpleParticles / 200);
    for (int particle = 0; particle < simpleParticles; particle++) {
        if (particle % report == 0)
            IJ.showProgress(particle, simpleParticles);
        double[] xyz = new double[3];
        if (linearDiffusion) {
            double[] dir = nextVector();
            for (int step = 0; step < simpleSteps; step++) {
                final double d = ((random[1].nextDouble() > 0.5) ? -1 : 1) * random[0].nextGaussian();
                for (int i = 0; i < 3; i++) {
                    xyz[i] += dir[i] * d;
                }
            }
        } else {
            for (int step = 0; step < simpleSteps; step++) {
                for (int i = 0; i < 3; i++) {
                    xyz[i] += random[i].nextGaussian();
                }
            }
        }
        for (int i = 0; i < 3; i++) xyz[i] *= scale;
        double msd = 0;
        for (int i = 0; i < 3; i++) {
            msd += xyz[i] * xyz[i];
            stats2[i].add(msd);
            // Store the actual distances
            stats[i].add(xyz[i]);
        }
    }
    IJ.showProgress(1);
    for (int i = 0; i < 3; i++) {
        plotJumpDistances(TITLE, stats2[i], i + 1);
        // Save stats to file for fitting
        save(stats2[i], i + 1, "msd");
        save(stats[i], i + 1, "d");
    }
    if (idCount > 0)
        new WindowOrganiser().tileWindows(idList);
}
Also used : StoredDataStatistics(gdsc.core.utils.StoredDataStatistics) WindowOrganiser(ij.plugin.WindowOrganiser) Well19937c(org.apache.commons.math3.random.Well19937c) RandomGenerator(org.apache.commons.math3.random.RandomGenerator)

Example 5 with WindowOrganiser

use of ij.plugin.WindowOrganiser in project GDSC-SMLM by aherbert.

the class TraceDiffusion method tileNewWindows.

private void tileNewWindows() {
    if (idCount > 0) {
        idList = Arrays.copyOf(idList, idCount);
        new WindowOrganiser().tileWindows(idList);
    }
}
Also used : WindowOrganiser(ij.plugin.WindowOrganiser)

Aggregations

WindowOrganiser (ij.plugin.WindowOrganiser)15 StoredDataStatistics (gdsc.core.utils.StoredDataStatistics)8 Statistics (gdsc.core.utils.Statistics)7 ArrayList (java.util.ArrayList)5 MemoryPeakResults (gdsc.smlm.results.MemoryPeakResults)3 ImageStack (ij.ImageStack)3 Plot2 (ij.gui.Plot2)3 PlotWindow (ij.gui.PlotWindow)3 Rectangle (java.awt.Rectangle)3 LinkedList (java.util.LinkedList)3 ArrayBlockingQueue (java.util.concurrent.ArrayBlockingQueue)3 ExecutorService (java.util.concurrent.ExecutorService)3 Future (java.util.concurrent.Future)3 BasePoint (gdsc.core.match.BasePoint)2 TurboList (gdsc.core.utils.TurboList)2 PeakResultPoint (gdsc.smlm.ij.plugins.ResultsMatchCalculator.PeakResultPoint)2 ImagePlus (ij.ImagePlus)2 ImageProcessor (ij.process.ImageProcessor)2 RandomGenerator (org.apache.commons.math3.random.RandomGenerator)2 Well19937c (org.apache.commons.math3.random.Well19937c)2