Search in sources :

Example 1 with AlignImagesFFT

use of gdsc.core.ij.AlignImagesFFT in project GDSC-SMLM by aherbert.

the class DriftCalculator method calculateUsingImageStack.

/**
	 * Calculates drift using images from a reference stack aligned to the overall z-projection.
	 * 
	 * @param stack
	 * 
	 * @param limits
	 * @return the drift { dx[], dy[] }
	 */
private double[][] calculateUsingImageStack(ImageStack stack, int[] limits) {
    // Update the limits using the stack size
    int upperT = startFrame + frameSpacing * (stack.getSize() - 1);
    limits[1] = FastMath.max(limits[1], upperT);
    // TODO - Truncate the stack if there are far too many frames for the localisation limits
    tracker.status("Constructing images");
    threadPool = Executors.newFixedThreadPool(Prefs.getThreads());
    // Built an image and FHT image for each slice
    final ImageProcessor[] images = new ImageProcessor[stack.getSize()];
    final FHT[] fhtImages = new FHT[stack.getSize()];
    List<Future<?>> futures = new LinkedList<Future<?>>();
    progressCounter = 0;
    totalCounter = images.length;
    int imagesPerThread = getImagesPerThread(images);
    final AlignImagesFFT aligner = new AlignImagesFFT();
    FloatProcessor referenceIp = stack.getProcessor(1).toFloat(0, null);
    // We do not care about the window method because this processor will not 
    // actually be used for alignment, it is a reference for the FHT size		
    aligner.init(referenceIp, WindowMethod.NONE, false);
    for (int i = 0; i < images.length; i += imagesPerThread) {
        futures.add(threadPool.submit(new ImageFHTInitialiser(stack, images, aligner, fhtImages, i, i + imagesPerThread)));
    }
    Utils.waitForCompletion(futures);
    tracker.progress(1);
    if (tracker.isEnded())
        return null;
    double[] dx = new double[limits[1] + 1];
    double[] dy = new double[dx.length];
    double[] originalDriftTimePoints = new double[dx.length];
    int[] blockT = new int[stack.getSize()];
    for (int i = 0, t = startFrame; i < stack.getSize(); i++, t += frameSpacing) {
        originalDriftTimePoints[t] = 1;
        blockT[i] = t;
    }
    double smoothing = updateSmoothingParameter(originalDriftTimePoints);
    lastdx = null;
    // For the first iteration calculate drift to the first image in the stack
    // (since the average projection may have a large drift blurring the image)
    double change = calculateDriftUsingImageStack(referenceIp, images, fhtImages, blockT, dx, dy, originalDriftTimePoints, smoothing, iterations);
    if (Double.isNaN(change) || tracker.isEnded())
        return null;
    plotDrift(limits, dx, dy);
    Utils.log("Drift Calculator : Initial drift " + Utils.rounded(change));
    for (int i = 1; i <= maxIterations; i++) {
        change = calculateDriftUsingImageStack(null, images, fhtImages, blockT, dx, dy, originalDriftTimePoints, smoothing, iterations);
        if (Double.isNaN(change))
            return null;
        plotDrift(limits, dx, dy);
        if (converged(i, change, getTotalDrift(dx, dy, originalDriftTimePoints)))
            break;
    }
    if (tracker.isEnded())
        return null;
    plotDrift(limits, dx, dy);
    return new double[][] { dx, dy };
}
Also used : FloatProcessor(ij.process.FloatProcessor) Point(java.awt.Point) LinkedList(java.util.LinkedList) ImageProcessor(ij.process.ImageProcessor) FHT(ij.process.FHT) Future(java.util.concurrent.Future) AlignImagesFFT(gdsc.core.ij.AlignImagesFFT)

Example 2 with AlignImagesFFT

use of gdsc.core.ij.AlignImagesFFT in project GDSC-SMLM by aherbert.

the class DriftCalculator method calculateDrift.

/**
	 * Calculate the drift of images to the reference image. Update the current drift parameters.
	 * 
	 * @param imageT
	 *            The frame number for each image
	 * @param scale
	 *            The image scale (used to adjust the drift to the correct size)
	 * @param dx
	 *            The X drift
	 * @param dy
	 *            The Y drift
	 * @param originalDriftTimePoints
	 *            Non-zero when the frame number refers to an aligned image frame
	 * @param smoothing
	 *            LOESS smoothing parameter
	 * @param iterations
	 *            LOESS iterations parameter
	 * @param images
	 *            The images to align
	 * @param reference
	 *            The reference image
	 * @param includeCurrentDrift
	 *            Set to true if the input images already have the current drift applied. The new drift will be added to
	 *            the current drift.
	 * @return
	 */
private double calculateDrift(int[] imageT, float scale, double[] dx, double[] dy, double[] originalDriftTimePoints, double smoothing, int iterations, final ImageProcessor[] images, FloatProcessor reference, boolean includeCurrentDrift) {
    // Align
    tracker.status("Aligning images");
    final AlignImagesFFT aligner = new AlignImagesFFT();
    aligner.init(reference, WindowMethod.NONE, false);
    final Rectangle alignBounds = AlignImagesFFT.createHalfMaxBounds(reference.getWidth(), reference.getHeight(), reference.getWidth(), reference.getHeight());
    List<double[]> alignments = Collections.synchronizedList(new ArrayList<double[]>(images.length));
    List<Future<?>> futures = new LinkedList<Future<?>>();
    int imagesPerThread = getImagesPerThread(images);
    for (int i = 0; i < images.length; i += imagesPerThread) {
        futures.add(threadPool.submit(new ImageAligner(aligner, images, imageT, alignBounds, alignments, i, i + imagesPerThread)));
    }
    Utils.waitForCompletion(futures);
    tracker.progress(1);
    // Used to flag when an alignment has failed
    originalDriftTimePoints = Arrays.copyOf(originalDriftTimePoints, originalDriftTimePoints.length);
    double[] newDx = new double[dx.length];
    double[] newDy = new double[dy.length];
    int ok = 0;
    for (double[] result : alignments) {
        int t = (int) result[2];
        if (Double.isNaN(result[0])) {
            // TODO: How to ignore bad alignments? 
            // Only do smoothing where there was an alignment?
            originalDriftTimePoints[t] = 0;
            tracker.log("WARNING : Unable to align image for time %d to the overall projection", t);
        } else {
            ok++;
            newDx[t] = result[0] / scale;
            newDy[t] = result[1] / scale;
            if (includeCurrentDrift) {
                // New drift = update + previous drift
                newDx[t] += dx[t];
                newDy[t] += dy[t];
            }
        }
    }
    if (ok < 2) {
        tracker.log("ERROR : Unable to align more than 1 image to the overall projection");
        return Double.NaN;
    }
    // Store the pure drift values for plotting
    calculatedTimepoints = Arrays.copyOf(originalDriftTimePoints, originalDriftTimePoints.length);
    lastdx = Arrays.copyOf(newDx, newDx.length);
    lastdy = Arrays.copyOf(newDy, newDy.length);
    // Perform smoothing
    if (smoothing > 0) {
        tracker.status("Smoothing drift");
        if (!smooth(newDx, newDy, originalDriftTimePoints, smoothing, iterations))
            return Double.NaN;
    }
    // Interpolate values for all time limits
    tracker.status("Interpolating drift");
    interpolate(newDx, newDy, originalDriftTimePoints);
    // Average drift correction for the calculated points should be zero to allow change comparison
    normalise(newDx, originalDriftTimePoints);
    normalise(newDy, originalDriftTimePoints);
    // Calculate change and update the input drift parameters
    double change = 0;
    for (int t = 0; t < dx.length; t++) {
        if (originalDriftTimePoints[t] != 0) {
            double d1 = dx[t] - newDx[t];
            double d2 = dy[t] - newDy[t];
            change += Math.sqrt(d1 * d1 + d2 * d2);
        }
        // Update all points since interpolation has already been done
        dx[t] = newDx[t];
        dy[t] = newDy[t];
    }
    tracker.status("");
    return change;
}
Also used : Rectangle(java.awt.Rectangle) AlignImagesFFT(gdsc.core.ij.AlignImagesFFT) Future(java.util.concurrent.Future) LinkedList(java.util.LinkedList) Point(java.awt.Point)

Aggregations

AlignImagesFFT (gdsc.core.ij.AlignImagesFFT)2 Point (java.awt.Point)2 LinkedList (java.util.LinkedList)2 Future (java.util.concurrent.Future)2 FHT (ij.process.FHT)1 FloatProcessor (ij.process.FloatProcessor)1 ImageProcessor (ij.process.ImageProcessor)1 Rectangle (java.awt.Rectangle)1