Examples with MemoryPeakResults uk.ac.sussex.gdsc.smlm.results.MemoryPeakResults used on opensource projects

Example 96 with MemoryPeakResults

use of uk.ac.sussex.gdsc.smlm.results.MemoryPeakResults in project GDSC-SMLM by aherbert.

the class BlinkEstimator method run.

@Override
public void run(String arg) {
    SmlmUsageTracker.recordPlugin(this.getClass(), arg);
    // Require some fit results and selected regions
    if (MemoryPeakResults.isMemoryEmpty()) {
        IJ.error(TITLE, "There are no fitting results in memory");
        return;
    }
    if (!showDialog()) {
        return;
    }
    final MemoryPeakResults results = ResultsManager.loadInputResults(settings.inputOption, true, null, null);
    if (MemoryPeakResults.isEmpty(results)) {
        IJ.error(TITLE, "No results could be loaded");
        IJ.showStatus("");
        return;
    }
    msPerFrame = results.getCalibrationReader().getExposureTime();
    ImageJUtils.log("%s: %d localisations", TITLE, results.size());
    showPlots = true;
    if (settings.rangeFittedPoints > 0) {
        computeFitCurves(results, true);
    } else {
        computeBlinkingRate(results, true);
    }
}

Example 97 with MemoryPeakResults

use of uk.ac.sussex.gdsc.smlm.results.MemoryPeakResults in project GDSC-SMLM by aherbert.

the class DensityEstimator method run.

@Override
public void run(String arg) {
    SmlmUsageTracker.recordPlugin(this.getClass(), arg);
    // Require some fit results and selected regions
    if (MemoryPeakResults.countMemorySize() == 0) {
        IJ.error(TITLE, "There are no fitting results in memory");
        return;
    }
    if (!showDialog()) {
        return;
    }
    // Currently this only supports pixel distance units
    final MemoryPeakResults results = ResultsManager.loadInputResults(settings.inputOption, false, DistanceUnit.PIXEL, null);
    if (MemoryPeakResults.isEmpty(results)) {
        IJ.error(TITLE, "No results could be loaded");
        IJ.showStatus("");
        return;
    }
    final long start = System.currentTimeMillis();
    IJ.showStatus("Calculating density ...");
    // Scale to um^2 from px^2
    final double scale = Math.pow(results.getDistanceConverter(DistanceUnit.UM).convertBack(1), 2);
    results.sort();
    final FrameCounter counter = results.newFrameCounter();
    final double localisationsPerFrame = (double) results.size() / (results.getLastFrame() - counter.currentFrame() + 1);
    final Rectangle bounds = results.getBounds(true);
    final double globalDensity = localisationsPerFrame / bounds.width / bounds.height;
    final int border = settings.border;
    final boolean includeSingles = settings.includeSingles;
    final int size = 2 * border + 1;
    final double minDensity = Math.pow(size, -2);
    ImageJUtils.log("%s : %s : Global density %s. Minimum density in %dx%d px = %s um^-2", TITLE, results.getName(), MathUtils.rounded(globalDensity * scale), size, size, MathUtils.rounded(minDensity * scale));
    final TIntArrayList x = new TIntArrayList();
    final TIntArrayList y = new TIntArrayList();
    final ExecutorService es = Executors.newFixedThreadPool(Prefs.getThreads());
    final LocalList<FrameDensity> densities = new LocalList<>();
    final LocalList<Future<?>> futures = new LocalList<>();
    results.forEach((PeakResultProcedure) (peak) -> {
        if (counter.advance(peak.getFrame())) {
            final FrameDensity fd = new FrameDensity(peak.getFrame(), x.toArray(), y.toArray(), border, includeSingles);
            densities.add(fd);
            futures.add(es.submit(fd));
            x.resetQuick();
            y.resetQuick();
        }
        x.add((int) peak.getXPosition());
        y.add((int) peak.getYPosition());
    });
    densities.add(new FrameDensity(counter.currentFrame(), x.toArray(), y.toArray(), border, includeSingles));
    futures.add(es.submit(densities.get(densities.size() - 1)));
    es.shutdown();
    // Wait
    ConcurrencyUtils.waitForCompletionUnchecked(futures);
    densities.sort((o1, o2) -> Integer.compare(o1.frame, o2.frame));
    final int total = densities.stream().mapToInt(fd -> fd.counts.length).sum();
    // Plot density
    final Statistics stats = new Statistics();
    final float[] frame = new float[total];
    final float[] density = new float[total];
    densities.stream().forEach(fd -> {
        for (int i = 0; i < fd.counts.length; i++) {
            final double d = (fd.counts[i] / fd.values[i]) * scale;
            frame[stats.getN()] = fd.frame;
            density[stats.getN()] = (float) d;
            stats.add(d);
        }
    });
    final double mean = stats.getMean();
    final double sd = stats.getStandardDeviation();
    final String label = String.format("Density = %s +/- %s um^-2", MathUtils.rounded(mean), MathUtils.rounded(sd));
    final Plot plot = new Plot("Frame vs Density", "Frame", "Density (um^-2)");
    plot.addPoints(frame, density, Plot.CIRCLE);
    plot.addLabel(0, 0, label);
    final WindowOrganiser wo = new WindowOrganiser();
    ImageJUtils.display(plot.getTitle(), plot, wo);
    // Histogram density
    new HistogramPlotBuilder("Local", StoredData.create(density), "Density (um^-2)").setPlotLabel(label).show(wo);
    wo.tile();
    // Log the number of singles
    final int singles = densities.stream().mapToInt(fd -> fd.singles).sum();
    ImageJUtils.log("Singles %d / %d (%s%%)", singles, results.size(), MathUtils.rounded(100.0 * singles / results.size()));
    IJ.showStatus(TITLE + " complete : " + TextUtils.millisToString(System.currentTimeMillis() - start));
}

Example 98 with MemoryPeakResults

use of uk.ac.sussex.gdsc.smlm.results.MemoryPeakResults in project GDSC-SMLM by aherbert.

the class DensityImage method cropWithBorder.

/**
 * Crop the results to the ROI. Add a border using the sampling radius so that counts do not have
 * to be approximated (i.e. all spots around the edge of the ROI will have a complete image to
 * sample from). The results are modified in place.
 *
 * @param results the results
 * @param isWithin Set to true if the added border is within the original bounds (i.e. no
 *        adjustment for missing counts is required)
 * @return the cropped results
 */
private MemoryPeakResults cropWithBorder(MemoryPeakResults results, boolean[] isWithin) {
    isWithin[0] = false;
    if (roiBounds == null) {
        return results;
    }
    // Adjust bounds relative to input results image:
    // Use the ROI relative to the frame the ROI is drawn on.
    // Map those fractional coordinates back to the original data bounds.
    final Rectangle bounds = results.getBounds();
    final double xscale = (double) roiImageWidth / bounds.width;
    final double yscale = (double) roiImageHeight / bounds.height;
    // Compute relative to the results bounds (if present)
    scaledRoiMinX = bounds.x + roiBounds.x / xscale;
    scaledRoiMaxX = scaledRoiMinX + roiBounds.width / xscale;
    scaledRoiMinY = bounds.y + roiBounds.y / yscale;
    scaledRoiMaxY = scaledRoiMinY + roiBounds.height / yscale;
    // Allow for the border
    final float minX = (int) (scaledRoiMinX - settings.radius);
    final float maxX = (int) Math.ceil(scaledRoiMaxX + settings.radius);
    final float minY = (int) (scaledRoiMinY - settings.radius);
    final float maxY = (int) Math.ceil(scaledRoiMaxY + settings.radius);
    // Create a new set of results within the bounds
    final MemoryPeakResults newResults = new MemoryPeakResults();
    newResults.begin();
    results.forEach(DistanceUnit.PIXEL, (XyrResultProcedure) (x, y, result) -> {
        if (x >= minX && x <= maxX && y >= minY && y <= maxY) {
            newResults.add(result);
        }
    });
    newResults.end();
    newResults.copySettings(results);
    newResults.setBounds(new Rectangle((int) minX, (int) minY, (int) (maxX - minX), (int) (maxY - minY)));
    isWithin[0] = minX >= bounds.x && minY >= bounds.y && maxX <= (bounds.x + bounds.width) && maxY <= (bounds.y + bounds.height);
    return newResults;
}

Example 99 with MemoryPeakResults

use of uk.ac.sussex.gdsc.smlm.results.MemoryPeakResults in project GDSC-SMLM by aherbert.

the class DiffusionRateTest method run.

@Override
public void run(String arg) {
    SmlmUsageTracker.recordPlugin(this.getClass(), arg);
    pluginSettings = Settings.load();
    pluginSettings.save();
    if (IJ.controlKeyDown()) {
        simpleTest();
        return;
    }
    extraOptions = ImageJUtils.isExtraOptions();
    if (!showDialog()) {
        return;
    }
    lastSimulation.set(null);
    final int totalSteps = (int) Math.ceil(settings.getSeconds() * settings.getStepsPerSecond());
    conversionFactor = 1000000.0 / (settings.getPixelPitch() * settings.getPixelPitch());
    // Diffusion rate is um^2/sec. Convert to pixels per simulation frame.
    final double diffusionRateInPixelsPerSecond = settings.getDiffusionRate() * conversionFactor;
    final double diffusionRateInPixelsPerStep = diffusionRateInPixelsPerSecond / settings.getStepsPerSecond();
    final double precisionInPixels = myPrecision / settings.getPixelPitch();
    final boolean addError = myPrecision != 0;
    ImageJUtils.log(TITLE + " : D = %s um^2/sec, Precision = %s nm", MathUtils.rounded(settings.getDiffusionRate(), 4), MathUtils.rounded(myPrecision, 4));
    ImageJUtils.log("Mean-displacement per dimension = %s nm/sec", MathUtils.rounded(1e3 * ImageModel.getRandomMoveDistance(settings.getDiffusionRate()), 4));
    if (extraOptions) {
        ImageJUtils.log("Step size = %s, precision = %s", MathUtils.rounded(ImageModel.getRandomMoveDistance(diffusionRateInPixelsPerStep)), MathUtils.rounded(precisionInPixels));
    }
    // Convert diffusion co-efficient into the standard deviation for the random walk
    final DiffusionType diffusionType = CreateDataSettingsHelper.getDiffusionType(settings.getDiffusionType());
    final double diffusionSigma = ImageModel.getRandomMoveDistance(diffusionRateInPixelsPerStep);
    ImageJUtils.log("Simulation step-size = %s nm", MathUtils.rounded(settings.getPixelPitch() * diffusionSigma, 4));
    // Move the molecules and get the diffusion rate
    IJ.showStatus("Simulating ...");
    final long start = System.nanoTime();
    final UniformRandomProvider random = UniformRandomProviders.create();
    final Statistics[] stats2D = new Statistics[totalSteps];
    final Statistics[] stats3D = new Statistics[totalSteps];
    final StoredDataStatistics jumpDistances2D = new StoredDataStatistics(totalSteps);
    final StoredDataStatistics jumpDistances3D = new StoredDataStatistics(totalSteps);
    for (int j = 0; j < totalSteps; j++) {
        stats2D[j] = new Statistics();
        stats3D[j] = new Statistics();
    }
    final SphericalDistribution dist = new SphericalDistribution(settings.getConfinementRadius() / settings.getPixelPitch());
    final Statistics asymptote = new Statistics();
    // Save results to memory
    final MemoryPeakResults results = new MemoryPeakResults(totalSteps);
    results.setCalibration(CalibrationHelper.create(settings.getPixelPitch(), 1, 1000.0 / settings.getStepsPerSecond()));
    results.setName(TITLE);
    results.setPsf(PsfHelper.create(PSFType.CUSTOM));
    int peak = 0;
    // Store raw coordinates
    final ArrayList<Point> points = new ArrayList<>(totalSteps);
    final StoredData totalJumpDistances1D = new StoredData(settings.getParticles());
    final StoredData totalJumpDistances2D = new StoredData(settings.getParticles());
    final StoredData totalJumpDistances3D = new StoredData(settings.getParticles());
    final NormalizedGaussianSampler gauss = SamplerUtils.createNormalizedGaussianSampler(random);
    for (int i = 0; i < settings.getParticles(); i++) {
        if (i % 16 == 0) {
            IJ.showProgress(i, settings.getParticles());
            if (ImageJUtils.isInterrupted()) {
                return;
            }
        }
        // Increment the frame so that tracing analysis can distinguish traces
        peak++;
        double[] origin = new double[3];
        final int id = i + 1;
        final MoleculeModel m = new MoleculeModel(id, origin.clone());
        if (addError) {
            origin = addError(origin, precisionInPixels, gauss);
        }
        if (pluginSettings.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 = (diffusionType == DiffusionType.LINEAR_WALK) ? nextVector(gauss) : null;
            for (int j = 0; j < totalSteps; j++) {
                double[] xyz = m.getCoordinates();
                final double[] originalXyz = xyz.clone();
                for (int n = pluginSettings.confinementAttempts; n-- > 0; ) {
                    if (diffusionType == DiffusionType.GRID_WALK) {
                        m.walk(diffusionSigma, random);
                    } else if (diffusionType == DiffusionType.LINEAR_WALK) {
                        m.slide(diffusionSigma, axis, random);
                    } 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, gauss);
                }
                peak = record(xyz, id, peak, stats2D[j], stats3D[j], jumpDistances2D, jumpDistances3D, origin, results);
            }
            asymptote.add(distance(m.getCoordinates()));
        } else if (diffusionType == 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, gauss);
                }
                peak = record(xyz, id, peak, stats2D[j], stats3D[j], jumpDistances2D, jumpDistances3D, origin, results);
            }
        } else if (diffusionType == DiffusionType.LINEAR_WALK) {
            final double[] axis = nextVector(gauss);
            for (int j = 0; j < totalSteps; j++) {
                m.slide(diffusionSigma, axis, random);
                double[] xyz = m.getCoordinates();
                points.add(new Point(id, xyz));
                if (addError) {
                    xyz = addError(xyz, precisionInPixels, gauss);
                }
                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, gauss);
                }
                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 long nanoseconds = System.nanoTime() - start;
    IJ.showProgress(1);
    MemoryPeakResults.addResults(results);
    simulation = new SimulationData(results.getName(), myPrecision);
    // Convert pixels^2/step to um^2/sec
    final double msd2D = (jumpDistances2D.getMean() / conversionFactor) / (results.getCalibrationReader().getExposureTime() / 1000);
    final double msd3D = (jumpDistances3D.getMean() / conversionFactor) / (results.getCalibrationReader().getExposureTime() / 1000);
    ImageJUtils.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", MathUtils.rounded(settings.getDiffusionRate()), MathUtils.rounded(myPrecision), jumpDistances2D.getN(), MathUtils.rounded(results.getCalibrationReader().getExposureTime() / 1000), MathUtils.rounded(jumpDistances2D.getMean() / conversionFactor), MathUtils.rounded(msd2D), MathUtils.rounded(jumpDistances3D.getMean() / conversionFactor), MathUtils.rounded(msd3D));
    aggregateIntoFrames(points, addError, precisionInPixels, gauss);
    IJ.showStatus("Analysing results ...");
    if (pluginSettings.showDiffusionExample) {
        showExample(totalSteps, diffusionSigma, random);
    }
    // Plot a graph of mean squared distance
    final double[] xValues = new double[stats2D.length];
    final double[] yValues2D = new double[stats2D.length];
    final double[] yValues3D = new double[stats3D.length];
    final double[] upper2D = new double[stats2D.length];
    final double[] lower2D = new double[stats2D.length];
    final double[] upper3D = new double[stats3D.length];
    final double[] lower3D = new double[stats3D.length];
    final SimpleRegression r2D = new SimpleRegression(false);
    final SimpleRegression r3D = new SimpleRegression(false);
    final int firstN = (pluginSettings.useConfinement) ? pluginSettings.fitN : totalSteps;
    for (int j = 0; j < totalSteps; j++) {
        // Convert steps to seconds
        xValues[j] = (j + 1) / settings.getStepsPerSecond();
        // 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;
    final PolynomialFunction 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 diffCoeff = best2D[1] / 4.0;
        final String msg = "2D Diffusion rate = " + MathUtils.rounded(diffCoeff, 4) + " um^2 / sec (" + TextUtils.nanosToString(nanoseconds) + ")";
        IJ.showStatus(msg);
        ImageJUtils.log(msg);
        diffCoeff = best3D[1] / 6.0;
        ImageJUtils.log("3D Diffusion rate = " + MathUtils.rounded(diffCoeff, 4) + " um^2 / sec (" + TextUtils.nanosToString(nanoseconds) + ")");
    } 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);
    // Show the total jump length for debugging
    // plotJumpDistances(TITLE + " total", totalJumpDistances1D, 1, totalSteps);
    // plotJumpDistances(TITLE + " total", totalJumpDistances2D, 2, totalSteps);
    // plotJumpDistances(TITLE + " total", totalJumpDistances3D, 3, totalSteps);
    windowOrganiser.tile();
    if (pluginSettings.useConfinement) {
        ImageJUtils.log("3D asymptote distance = %s nm (expected %.2f)", MathUtils.rounded(asymptote.getMean() * settings.getPixelPitch(), 4), 3 * settings.getConfinementRadius() / 4);
    }
}

Example 100 with MemoryPeakResults

use of uk.ac.sussex.gdsc.smlm.results.MemoryPeakResults in project GDSC-SMLM by aherbert.

the class CalibrateResults method showDialog.

private boolean showDialog(MemoryPeakResults results) {
    final ExtendedGenericDialog gd = new ExtendedGenericDialog(TITLE);
    gd.addHelp(HelpUrls.getUrl("calibrate-results"));
    final Calibration oldCalibration = results.getCalibration();
    final boolean existingCalibration = oldCalibration != null;
    if (!existingCalibration) {
        gd.addMessage("No calibration found, using defaults");
    }
    final CalibrationWriter cw = results.getCalibrationWriterSafe();
    gd.addStringField("Name", results.getName(), Math.max(Math.min(results.getName().length(), 60), 20));
    if (existingCalibration) {
        gd.addCheckbox("Update_all_linked_results", settings.updateAll);
    }
    PeakFit.addCameraOptions(gd, PeakFit.FLAG_QUANTUM_EFFICIENCY, cw);
    gd.addNumericField("Calibration (nm/px)", cw.getNmPerPixel(), 2);
    gd.addNumericField("Exposure_time (ms)", cw.getExposureTime(), 2);
    gd.showDialog();
    if (gd.wasCanceled()) {
        return false;
    }
    final String name = gd.getNextString();
    if (!results.getName().equals(name)) {
        // If the name is changed then remove and add back to memory
        final MemoryPeakResults existingResults = MemoryPeakResults.removeResults(results.getName());
        if (existingResults != null) {
            results = existingResults;
            results.setName(name);
            MemoryPeakResults.addResults(results);
        }
    }
    if (existingCalibration) {
        settings.updateAll = gd.getNextBoolean();
    }
    cw.setCameraType(SettingsManager.getCameraTypeValues()[gd.getNextChoiceIndex()]);
    cw.setNmPerPixel(Math.abs(gd.getNextNumber()));
    cw.setExposureTime(Math.abs(gd.getNextNumber()));
    gd.collectOptions();
    final Calibration newCalibration = cw.getCalibration();
    results.setCalibration(newCalibration);
    if (settings.updateAll) {
        // be missed.
        for (final MemoryPeakResults r : MemoryPeakResults.getAllResults()) {
            if (r.getCalibration() == oldCalibration) {
                r.setCalibration(newCalibration);
            }
        }
    }
    return true;
}