Search in sources :

Example 16 with Ticker

use of uk.ac.sussex.gdsc.core.logging.Ticker in project GDSC-SMLM by aherbert.

the class FailCountManager method analyseData.

private void analyseData() {
    final LocalList<FailCountData> failCountData = failCountDataRef.get();
    if (failCountData.isEmpty()) {
        IJ.error(TITLE, "No fail count data in memory");
        return;
    }
    if (!showAnalysisDialog()) {
        return;
    }
    final int maxCons = getMaxConsecutiveFailCount(failCountData);
    final int maxFail = getMaxFailCount(failCountData);
    // Create a set of fail counters
    final LocalList<FailCounter> counters = new LocalList<>();
    final TByteArrayList type = new TByteArrayList();
    for (int i = 0; i <= maxCons; i++) {
        counters.add(ConsecutiveFailCounter.create(i));
    }
    fill(type, counters, 0);
    // The other counters are user configured.
    // Ideally this would be a search to optimise the best parameters
    // for each counter as any enumeration may be way off the mark.
    // Note that 0 failures in a window can be scored using the consecutive fail counter.
    int max = Math.min(maxFail, settings.getRollingCounterMaxAllowedFailures());
    for (int fail = MathUtils.clip(1, maxFail, settings.getRollingCounterMinAllowedFailures()); fail <= max; fail++) {
        // Note that n-1 failures in window n can be scored using the consecutive fail counter.
        for (int window = Math.max(fail + 2, settings.getRollingCounterMinWindow()); window <= settings.getRollingCounterMaxWindow(); window++) {
            counters.add(RollingWindowFailCounter.create(fail, window));
        }
        switch(checkCounters(counters)) {
            case ANALYSE:
                break;
            case CONTINUE:
                break;
            case RETURN:
                return;
            default:
                throw new IllegalStateException();
        }
    }
    fill(type, counters, 1);
    max = Math.min(maxFail, settings.getWeightedCounterMaxAllowedFailures());
    for (int fail = MathUtils.min(maxFail, settings.getWeightedCounterMinAllowedFailures()); fail <= max; fail++) {
        for (int w = settings.getWeightedCounterMinPassDecrement(); w <= settings.getWeightedCounterMaxPassDecrement(); w++) {
            counters.add(WeightedFailCounter.create(fail, 1, w));
        }
        switch(checkCounters(counters)) {
            case ANALYSE:
                break;
            case CONTINUE:
                break;
            case RETURN:
                return;
            default:
                throw new IllegalStateException();
        }
    }
    fill(type, counters, 2);
    max = Math.min(maxFail, settings.getResettingCounterMaxAllowedFailures());
    for (int fail = MathUtils.min(maxFail, settings.getResettingCounterMinAllowedFailures()); fail <= max; fail++) {
        for (double f = settings.getResettingCounterMinResetFraction(); f <= settings.getResettingCounterMaxResetFraction(); f += settings.getResettingCounterIncResetFraction()) {
            counters.add(ResettingFailCounter.create(fail, f));
        }
        switch(checkCounters(counters)) {
            case ANALYSE:
                break;
            case CONTINUE:
                break;
            case RETURN:
                return;
            default:
                throw new IllegalStateException();
        }
    }
    fill(type, counters, 3);
    for (int count = settings.getPassRateCounterMinAllowedCounts(); count <= settings.getPassRateCounterMaxAllowedCounts(); count++) {
        for (double f = settings.getPassRateCounterMinPassRate(); f <= settings.getPassRateCounterMaxPassRate(); f += settings.getPassRateCounterIncPassRate()) {
            counters.add(PassRateFailCounter.create(count, f));
        }
        switch(checkCounters(counters)) {
            case ANALYSE:
                break;
            case CONTINUE:
                break;
            case RETURN:
                return;
            default:
                throw new IllegalStateException();
        }
    }
    fill(type, counters, 4);
    counters.trimToSize();
    // Score each of a set of standard fail counters against each frame using how
    // close they are to the target.
    final double[] score = new double[counters.size()];
    final double targetPassFraction = settings.getTargetPassFraction();
    final int nThreads = Prefs.getThreads();
    final ExecutorService executor = Executors.newFixedThreadPool(nThreads);
    final LocalList<Future<?>> futures = new LocalList<>(nThreads);
    final Ticker ticker = ImageJUtils.createTicker(failCountData.size(), nThreads);
    IJ.showStatus("Analysing " + TextUtils.pleural(counters.size(), "counter"));
    for (int i = 0; i < failCountData.size(); i++) {
        final FailCountData data = failCountData.unsafeGet(i);
        futures.add(executor.submit(() -> {
            if (IJ.escapePressed()) {
                return;
            }
            // TODO - Ideally this plugin should be run on benchmark data with ground truth.
            // The target could be to ensure all the correct results are fit
            // and false positives are excluded from incrementing the pass counter.
            // This could be done by saving the results from a benchmarking scoring
            // plugin to memory as the current dataset.
            data.initialiseAnalysis(targetPassFraction);
            // Score in blocks and then do a synchronized write to the combined score
            final Thread t = Thread.currentThread();
            final double[] s = new double[8192];
            int index = 0;
            while (index < counters.size()) {
                if (t.isInterrupted()) {
                    break;
                }
                final int block = Math.min(8192, counters.size() - index);
                for (int j = 0; j < block; j++) {
                    final FailCounter counter = counters.unsafeGet(index + j).newCounter();
                    s[j] = data.score(counter);
                }
                // Write to the combined score
                synchronized (score) {
                    for (int j = 0; j < block; j++) {
                        score[index + j] += s[j];
                    }
                }
                index += block;
            }
            ticker.tick();
        }));
    }
    executor.shutdown();
    ConcurrencyUtils.waitForCompletionUnchecked(futures);
    IJ.showProgress(1);
    if (IJ.escapePressed()) {
        IJ.showStatus("");
        IJ.error(TITLE, "Cancelled analysis");
        return;
    }
    IJ.showStatus("Summarising results ...");
    // TODO - check if the top filter is at the bounds of the range
    final int minIndex = SimpleArrayUtils.findMinIndex(score);
    ImageJUtils.log(TITLE + " Analysis : Best counter = %s (Score = %f)", counters.unsafeGet(minIndex).getDescription(), score[minIndex]);
    // Show a table of results for the top N for each type
    final int topN = Math.min(settings.getTableTopN(), score.length);
    if (topN > 0) {
        final byte[] types = type.toArray();
        final byte maxType = types[types.length - 1];
        final TextWindow resultsWindow = createTable();
        for (byte b = 0; b <= maxType; b++) {
            int[] indices;
            // Use a heap to avoid a full sort
            final IntDoubleMinHeap heap = new IntDoubleMinHeap(topN);
            for (int i = 0; i < score.length; i++) {
                if (types[i] == b) {
                    heap.offer(score[i], i);
                }
            }
            if (heap.getSize() == 0) {
                continue;
            }
            indices = heap.getItems();
            // Ensure sorted
            SortUtils.sortIndices(indices, score, false);
            final StringBuilder sb = new StringBuilder();
            try (final BufferedTextWindow tw = new BufferedTextWindow(resultsWindow)) {
                for (int i = 0; i < topN; i++) {
                    sb.setLength(0);
                    final int j = indices[i];
                    sb.append(i + 1).append('\t');
                    sb.append(counters.unsafeGet(j).getDescription()).append('\t');
                    sb.append(score[j]);
                    tw.append(sb.toString());
                }
            }
        }
    }
    // TODO - Save the best fail counter to the current fit configuration.
    IJ.showStatus("");
}
Also used : RollingWindowFailCounter(uk.ac.sussex.gdsc.smlm.results.count.RollingWindowFailCounter) FailCounter(uk.ac.sussex.gdsc.smlm.results.count.FailCounter) ResettingFailCounter(uk.ac.sussex.gdsc.smlm.results.count.ResettingFailCounter) ConsecutiveFailCounter(uk.ac.sussex.gdsc.smlm.results.count.ConsecutiveFailCounter) WeightedFailCounter(uk.ac.sussex.gdsc.smlm.results.count.WeightedFailCounter) PassRateFailCounter(uk.ac.sussex.gdsc.smlm.results.count.PassRateFailCounter) TByteArrayList(gnu.trove.list.array.TByteArrayList) BufferedTextWindow(uk.ac.sussex.gdsc.core.ij.BufferedTextWindow) Ticker(uk.ac.sussex.gdsc.core.logging.Ticker) LocalList(uk.ac.sussex.gdsc.core.utils.LocalList) TextWindow(ij.text.TextWindow) BufferedTextWindow(uk.ac.sussex.gdsc.core.ij.BufferedTextWindow) IntDoubleMinHeap(uk.ac.sussex.gdsc.core.trees.heaps.IntDoubleMinHeap) ExecutorService(java.util.concurrent.ExecutorService) Future(java.util.concurrent.Future)

Example 17 with Ticker

use of uk.ac.sussex.gdsc.core.logging.Ticker in project GDSC-SMLM by aherbert.

the class CubicSplineManager method createCubicSpline.

/**
 * Creates the cubic spline.
 *
 * @param imagePsf the image PSF details
 * @param image the image
 * @param singlePrecision Set to true to use single precision (float values) to store the cubic
 *        spline coefficients
 * @return the cubic spline PSF
 */
public static CubicSplinePsf createCubicSpline(ImagePSFOrBuilder imagePsf, ImageStack image, final boolean singlePrecision) {
    final int maxx = image.getWidth();
    final int maxy = image.getHeight();
    final int maxz = image.getSize();
    final float[][] psf = new float[maxz][];
    for (int z = 0; z < maxz; z++) {
        psf[z] = ImageJImageConverter.getData(image.getPixels(z + 1), null);
    }
    // We reduce by a factor of 3
    final int maxi = (maxx - 1) / 3;
    final int maxj = (maxy - 1) / 3;
    final int maxk = (maxz - 1) / 3;
    final int size = maxi * maxj;
    final CustomTricubicFunction[][] splines = new CustomTricubicFunction[maxk][size];
    final int threadCount = Prefs.getThreads();
    final Ticker ticker = ImageJUtils.createTicker((long) maxi * maxj * maxk, threadCount);
    final ExecutorService threadPool = Executors.newFixedThreadPool(threadCount);
    final LocalList<Future<?>> futures = new LocalList<>(maxk);
    // spline node along each dimension, i.e. dimension length = n*3 + 1 with n the number of nodes.
    for (int k = 0; k < maxk; k++) {
        final int kk = k;
        futures.add(threadPool.submit(() -> {
            final CubicSplineCalculator calc = new CubicSplineCalculator();
            final double[] value = new double[64];
            final int zz = 3 * kk;
            for (int j = 0, index = 0; j < maxj; j++) {
                // 4x4 block origin in the XY data
                int index0 = 3 * j * maxx;
                for (int i = 0; i < maxi; i++, index++) {
                    ticker.tick();
                    int count = 0;
                    for (int z = 0; z < 4; z++) {
                        final float[] data = psf[zz + z];
                        for (int y = 0; y < 4; y++) {
                            for (int x = 0, ii = index0 + y * maxx; x < 4; x++) {
                                value[count++] = data[ii++];
                            }
                        }
                    }
                    splines[kk][index] = CustomTricubicFunctionUtils.create(calc.compute(value));
                    if (singlePrecision) {
                        splines[kk][index] = splines[kk][index].toSinglePrecision();
                    }
                    index0 += 3;
                }
            }
        }));
    }
    ticker.stop();
    threadPool.shutdown();
    ConcurrencyUtils.waitForCompletionUnchecked(futures);
    // Normalise
    double maxSum = 0;
    for (int k = 0; k < maxk; k++) {
        double sum = 0;
        for (int i = 0; i < size; i++) {
            sum += splines[k][i].value000();
        }
        if (maxSum < sum) {
            maxSum = sum;
        }
    }
    if (maxSum == 0) {
        throw new IllegalStateException("The cubic spline has no maximum signal");
    }
    final double scale = 1.0 / maxSum;
    for (int k = 0; k < maxk; k++) {
        for (int i = 0; i < size; i++) {
            splines[k][i] = splines[k][i].scale(scale);
        }
    }
    // Create on an integer scale
    final CubicSplineData f = new CubicSplineData(maxi, maxj, splines);
    // Create a new info with the PSF details
    final ImagePSF.Builder b = ImagePSF.newBuilder();
    b.setImageCount(imagePsf.getImageCount());
    // Reducing the image has the effect of enlarging the pixel size
    b.setPixelSize(imagePsf.getPixelSize() * 3.0);
    b.setPixelDepth(imagePsf.getPixelDepth() * 3.0);
    // The centre has to be moved as we reduced the image size by 3.
    // In the ImagePSF the XY centre puts 0.5 at the centre of the pixel.
    // The spline puts 0,0 at the centre of each pixel for convenience.
    double cx = maxi / 2.0;
    if (imagePsf.getXCentre() != 0) {
        cx = (imagePsf.getXCentre() - 0.5) / 3;
    }
    double cy = maxj / 2.0;
    if (imagePsf.getYCentre() != 0) {
        cy = (imagePsf.getYCentre() - 0.5) / 3;
    }
    double cz = maxk / 2.0;
    if (imagePsf.getZCentre() != 0) {
        cz = imagePsf.getZCentre() / 3;
    } else if (imagePsf.getCentreImage() != 0) {
        cz = (imagePsf.getCentreImage() - 1) / 3.0;
    }
    b.setXCentre(cx);
    b.setYCentre(cy);
    b.setZCentre(cz);
    return new CubicSplinePsf(b.build(), f);
}
Also used : CubicSplineCalculator(uk.ac.sussex.gdsc.smlm.function.cspline.CubicSplineCalculator) Ticker(uk.ac.sussex.gdsc.core.logging.Ticker) LocalList(uk.ac.sussex.gdsc.core.utils.LocalList) ImagePSF(uk.ac.sussex.gdsc.smlm.data.config.PSFProtos.ImagePSF) CubicSplineData(uk.ac.sussex.gdsc.smlm.function.cspline.CubicSplineData) ExecutorService(java.util.concurrent.ExecutorService) Future(java.util.concurrent.Future) CustomTricubicFunction(uk.ac.sussex.gdsc.core.math.interpolation.CustomTricubicFunction)

Example 18 with Ticker

use of uk.ac.sussex.gdsc.core.logging.Ticker in project GDSC-SMLM by aherbert.

the class AstigmatismModelManager method fitRegion.

private boolean fitRegion() {
    final int radius = config.getFittingWidth();
    if (pluginSettings.getLogFitProgress()) {
        fitConfig.setLog(ImageJPluginLoggerHelper.getLogger(getClass()));
    }
    // Create a fit engine
    results = new MemoryPeakResults();
    results.setCalibration(fitConfig.getCalibration());
    results.setPsf(fitConfig.getPsf());
    results.setSortAfterEnd(true);
    results.begin();
    final int threadCount = Prefs.getThreads();
    final FitEngine engine = FitEngine.create(config, SynchronizedPeakResults.create(results, threadCount), threadCount, FitQueue.BLOCKING);
    final IJImageSource source = new IJImageSource(imp);
    source.open();
    final Rectangle r1 = new Rectangle(cx - radius, cy - radius, 2 * radius + 1, 2 * radius + 1);
    final Rectangle regionBounds = r1.intersection(new Rectangle(source.getWidth(), source.getHeight()));
    // Fit only a spot in the centre
    final int x = cx - regionBounds.x;
    final int y = cy - regionBounds.y;
    final int[] maxIndices = new int[] { y * regionBounds.width + x };
    final Ticker ticker = ImageJUtils.createTicker(source.getFrames(), threadCount);
    IJ.showStatus("Fitting ...");
    boolean shutdown = false;
    while (!shutdown) {
        // Extract the region from each frame
        final float[] region = source.next(regionBounds);
        if (region == null) {
            break;
        }
        final FitParameters params = new FitParameters();
        params.maxIndices = maxIndices.clone();
        final int slice = (int) ticker.getCurrent();
        final ParameterisedFitJob job = new ParameterisedFitJob(slice, params, slice, region, regionBounds);
        engine.run(job);
        ticker.tick();
        shutdown = IJ.escapePressed();
    }
    if (shutdown) {
        IJ.showStatus("Cancelled");
    }
    engine.end(shutdown);
    results.end();
    IJ.showProgress(1);
    if (!shutdown) {
        ImageJUtils.log("Fit %d/%s", results.size(), TextUtils.pleural(source.getFrames(), "spot"));
    }
    return !shutdown;
}
Also used : FitParameters(uk.ac.sussex.gdsc.smlm.engine.FitParameters) IJImageSource(uk.ac.sussex.gdsc.smlm.ij.IJImageSource) FitEngine(uk.ac.sussex.gdsc.smlm.engine.FitEngine) ParameterisedFitJob(uk.ac.sussex.gdsc.smlm.engine.ParameterisedFitJob) Ticker(uk.ac.sussex.gdsc.core.logging.Ticker) Rectangle(java.awt.Rectangle) MemoryPeakResults(uk.ac.sussex.gdsc.smlm.results.MemoryPeakResults)

Example 19 with Ticker

use of uk.ac.sussex.gdsc.core.logging.Ticker in project GDSC-SMLM by aherbert.

the class CubicSplineData method sample.

/**
 * Sample the function.
 *
 * <p>n samples will be taken per node in each dimension. A final sample is taken at the end of
 * the sample range thus the final range for each axis will be the current axis range.
 *
 * <p>The procedure setValue(int,int,int,double) method will be executed in ZYX order.
 *
 * @param nx the number of samples per spline node in the x dimension
 * @param ny the number of samples per spline node in the y dimension
 * @param nz the number of samples per spline node in the z dimension
 * @param procedure the procedure
 * @param progress the progress
 * @throws IllegalArgumentException If the number of sample is not positive
 */
public void sample(int nx, int ny, int nz, TrivalueProcedure procedure, TrackProgress progress) {
    if (nx < 1 || ny < 1 || nz < 1) {
        throw new IllegalArgumentException("Samples must be positive");
    }
    // We can interpolate all nodes n-times plus a final point at the last node
    final int maxx = (getMaxX() - 1) * nx;
    final int maxy = (getMaxY() - 1) * ny;
    final int maxz = (getMaxZ() - 1) * nz;
    if (!procedure.setDimensions(maxx + 1, maxy + 1, maxz + 1)) {
        return;
    }
    final Ticker ticker = Ticker.createStarted(progress, (long) (maxx + 1) * (maxy + 1) * (maxz + 1), false);
    // Pre-compute interpolation tables
    final CubicSplinePosition[] sx = createCubicSplinePosition(nx);
    final CubicSplinePosition[] sy = createCubicSplinePosition(ny);
    final CubicSplinePosition[] sz = createCubicSplinePosition(nz);
    // Write axis values
    // Cache the table and the spline position to use for each interpolation point
    final int[] xt = new int[maxx + 1];
    final int[] xp = new int[maxx + 1];
    for (int x = 0; x <= maxx; x++) {
        int xposition = x / nx;
        int xtable = x % nx;
        if (x == maxx) {
            // Final interpolation point
            xposition--;
            xtable = nx;
        }
        xt[x] = xtable;
        xp[x] = xposition;
        procedure.setX(x, xposition + (double) xtable / nx);
    }
    final int[] yt = new int[maxy + 1];
    final int[] yp = new int[maxy + 1];
    for (int y = 0; y <= maxy; y++) {
        int yposition = y / ny;
        int ytable = y % ny;
        if (y == maxy) {
            // Final interpolation point
            yposition--;
            ytable = ny;
        }
        yt[y] = ytable;
        yp[y] = yposition;
        procedure.setY(y, yposition + (double) ytable / ny);
    }
    final int[] zt = new int[maxz + 1];
    final int[] zp = new int[maxz + 1];
    for (int z = 0; z <= maxz; z++) {
        int zposition = z / nz;
        int ztable = z % nz;
        if (z == maxz) {
            // Final interpolation point
            zposition--;
            ztable = nz;
        }
        zt[z] = ztable;
        zp[z] = zposition;
        procedure.setZ(z, zposition + (double) ztable / nz);
    }
    // Write interpolated values
    for (int z = 0; z <= maxz; z++) {
        final CustomTricubicFunction[] xySplines = splines[zp[z]];
        for (int y = 0; y <= maxy; y++) {
            final int index = yp[y] * getMaxX();
            for (int x = 0; x <= maxx; x++) {
                procedure.setValue(x, y, z, xySplines[index + xp[x]].value(sx[xt[x]], sy[yt[y]], sz[zt[z]]));
                ticker.tick();
            }
        }
    }
    ticker.stop();
}
Also used : Ticker(uk.ac.sussex.gdsc.core.logging.Ticker) CubicSplinePosition(uk.ac.sussex.gdsc.core.math.interpolation.CubicSplinePosition) FloatCustomTricubicFunction(uk.ac.sussex.gdsc.core.math.interpolation.FloatCustomTricubicFunction) CustomTricubicFunction(uk.ac.sussex.gdsc.core.math.interpolation.CustomTricubicFunction)

Example 20 with Ticker

use of uk.ac.sussex.gdsc.core.logging.Ticker in project GDSC-SMLM by aherbert.

the class CubicSplineData method write.

/**
 * Write a tricubic splines to the output stream. The output will be buffered for performance.
 *
 * @param outputStream the output stream
 * @param progress the progress
 * @throws IOException Signals that an I/O exception has occurred.
 */
public void write(OutputStream outputStream, TrackProgress progress) throws IOException {
    // Write dimensions
    final int maxz = splines.length;
    final Ticker ticker = Ticker.createStarted(progress, (long) maxx * maxy * maxz, false);
    final BufferedOutputStream buffer = new BufferedOutputStream(outputStream);
    final DataOutput out = new DataOutputStream(buffer);
    out.writeInt(maxx);
    out.writeInt(maxy);
    out.writeInt(maxz);
    // Write precision
    final boolean singlePrecision = isSinglePrecision();
    out.writeBoolean(singlePrecision);
    final SplineWriter writer = (singlePrecision) ? new FloatSplineWriter() : new DoubleSplineWriter();
    final int size = maxx * maxy;
    for (int z = 0; z < maxz; z++) {
        for (int i = 0; i < size; i++) {
            writer.write(out, splines[z][i]);
            ticker.tick();
        }
    }
    buffer.flush();
    ticker.stop();
}
Also used : DataOutput(java.io.DataOutput) DataOutputStream(java.io.DataOutputStream) Ticker(uk.ac.sussex.gdsc.core.logging.Ticker) BufferedOutputStream(java.io.BufferedOutputStream)

Aggregations

Ticker (uk.ac.sussex.gdsc.core.logging.Ticker)33 LinkedList (java.util.LinkedList)16 Future (java.util.concurrent.Future)12 ImageStack (ij.ImageStack)11 ArrayBlockingQueue (java.util.concurrent.ArrayBlockingQueue)10 LocalList (uk.ac.sussex.gdsc.core.utils.LocalList)10 ConcurrentRuntimeException (org.apache.commons.lang3.concurrent.ConcurrentRuntimeException)9 ExecutorService (java.util.concurrent.ExecutorService)8 AtomicInteger (java.util.concurrent.atomic.AtomicInteger)8 ImagePlus (ij.ImagePlus)7 ArrayList (java.util.ArrayList)6 List (java.util.List)6 Statistics (uk.ac.sussex.gdsc.core.utils.Statistics)6 MemoryPeakResults (uk.ac.sussex.gdsc.smlm.results.MemoryPeakResults)6 TIntObjectHashMap (gnu.trove.map.hash.TIntObjectHashMap)5 ImageProcessor (ij.process.ImageProcessor)5 TextWindow (ij.text.TextWindow)5 Point (java.awt.Point)5 Rectangle (java.awt.Rectangle)5 UniformRandomProvider (org.apache.commons.rng.UniformRandomProvider)5