Examples with IntHistogram uk.ac.sussex.gdsc.core.threshold.IntHistogram used on opensource projects

Example 1 with IntHistogram

use of uk.ac.sussex.gdsc.core.threshold.IntHistogram in project GDSC-SMLM by aherbert.

the class CameraModelAnalysis method execute.

/**
 * Execute the analysis.
 */
private boolean execute() {
    dirty = false;
    final CameraModelAnalysisSettings settings = this.settings.build();
    if (!(getGain(settings) > 0)) {
        ImageJUtils.log(TITLE + "Error: No total gain");
        return false;
    }
    if (!(settings.getPhotons() > 0)) {
        ImageJUtils.log(TITLE + "Error: No photons");
        return false;
    }
    // Avoid repeating the same analysis
    if (settings.equals(lastSettings)) {
        return true;
    }
    lastSettings = settings;
    final IntHistogram h = getHistogram(settings);
    // Build cumulative distribution
    final double[][] cdf1 = cumulativeHistogram(h);
    final double[] x1 = cdf1[0];
    final double[] y1 = cdf1[1];
    // Interpolate to 300 steps faster evaluation?
    // Get likelihood function
    final LikelihoodFunction f = getLikelihoodFunction(settings);
    // Create likelihood cumulative distribution
    final double[][] cdf2 = cumulativeDistribution(settings, cdf1, f);
    // Compute Komolgorov distance
    final double[] distanceAndValue = getDistance(cdf1, cdf2);
    final double distance = distanceAndValue[0];
    final double value = distanceAndValue[1];
    final double area = distanceAndValue[2];
    final double[] x2 = cdf2[0];
    final double[] y2 = cdf2[1];
    // Fill y1
    int offset = 0;
    while (x2[offset] < x1[0]) {
        offset++;
    }
    final double[] y1b = new double[y2.length];
    System.arraycopy(y1, 0, y1b, offset, y1.length);
    Arrays.fill(y1b, offset + y1.length, y2.length, y1[y1.length - 1]);
    // KolmogorovSmirnovTest
    // n is the number of samples used to build the probability distribution.
    final int n = (int) MathUtils.sum(h.histogramCounts);
    // From KolmogorovSmirnovTest.kolmogorovSmirnovTest(RealDistribution distribution, double[]
    // data, boolean exact):
    // Returns the p-value associated with the null hypothesis that data is a sample from
    // distribution.
    // E.g. If p<0.05 then the null hypothesis is rejected and the data do not match the
    // distribution.
    double pvalue = Double.NaN;
    try {
        pvalue = 1d - kolmogorovSmirnovTest.cdf(distance, n);
    } catch (final MathArithmeticException ex) {
    // Cannot be computed to leave at NaN
    }
    // Plot
    final WindowOrganiser wo = new WindowOrganiser();
    String title = TITLE + " CDF";
    Plot plot = new Plot(title, "Count", "CDF");
    final double max = 1.05 * MathUtils.maxDefault(1, y2);
    plot.setLimits(x2[0], x2[x2.length - 1], 0, max);
    plot.setColor(Color.blue);
    plot.addPoints(x2, y1b, Plot.BAR);
    plot.setColor(Color.red);
    plot.addPoints(x2, y2, Plot.BAR);
    plot.setColor(Color.magenta);
    plot.drawLine(value, 0, value, max);
    plot.setColor(Color.black);
    plot.addLegend("CDF\nModel");
    plot.addLabel(0, 0, String.format("Distance=%s @ %.0f (Mean=%s) : p=%s", MathUtils.rounded(distance), value, MathUtils.rounded(area), MathUtils.rounded(pvalue)));
    ImageJUtils.display(title, plot, ImageJUtils.NO_TO_FRONT, wo);
    // Show the histogram
    title = TITLE + " Histogram";
    plot = new Plot(title, "Count", "Frequency");
    plot.setLimits(x1[0] - 0.5, x1[x1.length - 1] + 1.5, 0, MathUtils.max(h.histogramCounts) * 1.05);
    plot.setColor(Color.blue);
    plot.addPoints(x1, SimpleArrayUtils.toDouble(h.histogramCounts), Plot.BAR);
    plot.setColor(Color.red);
    final double[] x = floatHistogram[0].clone();
    final double[] y = floatHistogram[1].clone();
    final double scale = n / (MathUtils.sum(y) * (x[1] - x[0]));
    for (int i = 0; i < y.length; i++) {
        y[i] *= scale;
    }
    plot.addPoints(x, y, Plot.LINE);
    plot.setColor(Color.black);
    plot.addLegend("Sample\nExpected");
    ImageJUtils.display(title, plot, ImageJUtils.NO_TO_FRONT, wo);
    wo.tile();
    return true;
}