Examples with PeakResult gdsc.smlm.results.PeakResult used on opensource projects

Example 76 with PeakResult

use of gdsc.smlm.results.PeakResult in project GDSC-SMLM by aherbert.

the class FitWorker method findNeighbours.

/**
	 * Search for any neighbours within a set height of the specified peak that is within the search region bounds.
	 *
	 * @param regionBounds
	 *            the region bounds
	 * @param n
	 *            the candidate index
	 * @param background
	 *            The background in the region
	 * @return The number of neighbours
	 */
private int findNeighbours(Rectangle regionBounds, int n, float background) {
    final int xmin = regionBounds.x;
    final int xmax = xmin + regionBounds.width - 1;
    final int ymin = regionBounds.y;
    final int ymax = ymin + regionBounds.height - 1;
    final Candidate spot = candidates.get(n);
    final float heightThreshold;
    if (relativeIntensity) {
        // No background when spot filter has relative intensity
        heightThreshold = (float) (spot.intensity * config.getNeighbourHeightThreshold());
    } else {
        if (spot.intensity < background)
            heightThreshold = spot.intensity;
        else
            heightThreshold = (float) ((spot.intensity - background) * config.getNeighbourHeightThreshold() + background);
    }
    // Check all maxima that are lower than this
    candidateNeighbourCount = 0;
    // Using the neighbour grid.
    // Note this will also include all higher intensity spots that failed to be fit. 
    // These may still have estimates.
    final CandidateList neighbours = findNeighbours(spot);
    for (int i = 0; i < neighbours.size; i++) {
        final Candidate neighbour = neighbours.get(i);
        if (isFit(neighbour.index))
            continue;
        if (canIgnore(neighbour.x, neighbour.y, xmin, xmax, ymin, ymax, neighbour.intensity, heightThreshold))
            continue;
        candidateNeighbours[candidateNeighbourCount++] = neighbour;
    }
    // XXX Debugging
    //		int c = 0;
    //		// Processing all lower spots.
    //		//for (int i = n + 1; i < candidates.getSize(); i++)
    //		// Processing all spots.
    //		for (int i = 0; i < this.candidates.getSize(); i++)
    //		{
    //			if (i == n || isFit(i))
    //				continue;
    //			if (canIgnore(this.candidates.get(i).x, this.candidates.get(i).y, xmin, xmax, ymin, ymax,
    //					this.candidates.get(i).intensity, heightThreshold))
    //				continue;
    //			//neighbourIndices[c++] = i;
    //			if (neighbourIndices[c++] != i)
    //				throw new RuntimeException("invalid grid neighbours");
    //		}
    // Check all existing maxima. 
    fittedNeighbourCount = 0;
    if (!sliceResults.isEmpty()) {
        // Since these will be higher than the current peak it is prudent to extend the range that should be considered.
        // Use 2x the configured peak standard deviation.
        //final float range = 2f *
        //		(float) FastMath.max(fitConfig.getInitialPeakStdDev0(), fitConfig.getInitialPeakStdDev1());
        //final float xmin2 = regionBounds.x - range;
        //final float xmax2 = regionBounds.x + regionBounds.width + range;
        //final float ymin2 = regionBounds.y - range;
        //final float ymax2 = regionBounds.y + regionBounds.height + range;
        //for (int i = 0; i < sliceResults.size(); i++)
        //{
        //	final PeakResult result = sliceResults.get(i);
        //	// No height threshold check as this is a validated peak
        //	if (canIgnore(result.getXPosition(), result.getYPosition(), xmin2, xmax2, ymin2, ymax2))
        //		continue;
        //	fittedNeighbourIndices[fittedNeighbourCount++] = i;
        //}
        // Note: A smarter filter would be to compute the bounding rectangle of each fitted result and see if it 
        // overlaps the target region. This would involve overlap analysis
        final double x0min = regionBounds.x;
        final double y0min = regionBounds.y;
        final double x0max = regionBounds.x + regionBounds.width;
        final double y0max = regionBounds.y + regionBounds.height;
        final PeakResult[] peakNeighbours = findPeakNeighbours(spot);
        for (int i = 0; i < peakNeighbours.length; i++) {
            final PeakResult neighbour = peakNeighbours[i];
            // No height threshold check as this is a validated peak
            final double xw = 2 * neighbour.getXSD();
            final double yw = 2 * neighbour.getYSD();
            if (intersects(x0min, y0min, x0max, y0max, neighbour.getXPosition() - xw, neighbour.getYPosition() - yw, neighbour.getXPosition() + xw, neighbour.getYPosition() + yw))
                fittedNeighbours[fittedNeighbourCount++] = neighbour;
        }
    }
    return candidateNeighbourCount + fittedNeighbourCount;
}

Example 77 with PeakResult

use of gdsc.smlm.results.PeakResult in project GDSC-SMLM by aherbert.

the class OptimumDistanceResultFilter method filter.

/*
	 * (non-Javadoc)
	 * 
	 * @see gdsc.smlm.engine.filter.ResultFilter#filter(gdsc.smlm.fitting.FitResult, int,
	 * gdsc.smlm.results.PeakResult[])
	 */
@Override
public void filter(FitResult fitResult, int maxIndex, PeakResult... results) {
    for (PeakResult r : results) {
        if (r == null)
            continue;
        for (int i = 0; i < filter.size(); i++) {
            float[] coord = filter.get(i);
            final float dx = r.getXPosition() - coord[0];
            final float dy = r.getYPosition() - coord[1];
            // Only check if within the distance threshold 
            if (dx * dx + dy * dy < d2) {
                // Then filter by signal strength
                float s = r.getSignal();
                if (s < bestSignal[i])
                    continue;
                bestFitResults[i] = fitResult;
                bestIndices[i] = maxIndex;
                bestSignal[i] = s;
                bestPeakResults[i] = r;
            }
        }
    }
}

Example 78 with PeakResult

use of gdsc.smlm.results.PeakResult in project GDSC-SMLM by aherbert.

the class DriftCalculator method calculateUsingFrames.

/**
	 * Calculates drift using images from N consecutive frames aligned to the overall image.
	 * 
	 * @param results
	 * @param limits
	 * @param reconstructionSize
	 * @return the drift { dx[], dy[] }
	 */
private double[][] calculateUsingFrames(MemoryPeakResults results, int[] limits, int reconstructionSize) {
    double[] dx = new double[limits[1] + 1];
    double[] dy = new double[dx.length];
    // Extract the localisations into blocks of N consecutive frames
    ArrayList<ArrayList<Localisation>> blocks = new ArrayList<ArrayList<Localisation>>();
    results.sort();
    int t = 0;
    ArrayList<Localisation> nextBlock = null;
    for (PeakResult r : results.getResults()) {
        if (r.getFrame() > t) {
            while (r.getFrame() > t) t += frames;
            // To avoid blocks without many results only create a new block if the min size has been met
            if (nextBlock == null || nextBlock.size() >= minimimLocalisations)
                nextBlock = new ArrayList<Localisation>();
            blocks.add(nextBlock);
        }
        nextBlock.add(new Localisation(r.getFrame(), r.getXPosition(), r.getYPosition(), r.getSignal()));
    }
    if (blocks.size() < 2) {
        tracker.log("ERROR : Require at least 2 images for drift calculation");
        return null;
    }
    // Check the final block has enough localisations
    if (nextBlock.size() < minimimLocalisations) {
        blocks.remove(blocks.size() - 1);
        ArrayList<Localisation> combinedBlock = blocks.get(blocks.size() - 1);
        combinedBlock.addAll(nextBlock);
        if (blocks.size() < 2) {
            tracker.log("ERROR : Require at least 2 images for drift calculation");
            return null;
        }
    }
    // Find the average time point for each block
    int[] blockT = new int[blocks.size()];
    t = 0;
    for (ArrayList<Localisation> block : blocks) {
        long sum = 0;
        for (Localisation r : block) {
            sum += r.t;
        }
        blockT[t++] = (int) (sum / block.size());
    }
    // Calculate a scale to use when constructing the images for alignment
    Rectangle bounds = results.getBounds(true);
    float scale = (reconstructionSize - 1f) / FastMath.max(bounds.width, bounds.height);
    threadPool = Executors.newFixedThreadPool(Prefs.getThreads());
    double[] originalDriftTimePoints = getOriginalDriftTimePoints(dx, blockT);
    lastdx = null;
    double smoothing = updateSmoothingParameter(originalDriftTimePoints);
    double change = calculateDriftUsingFrames(blocks, blockT, bounds, scale, 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 = calculateDriftUsingFrames(blocks, blockT, bounds, scale, 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 };
}

Example 79 with PeakResult

use of gdsc.smlm.results.PeakResult in project GDSC-SMLM by aherbert.

the class DriftCalculator method findSpots.

private Spot[] findSpots(MemoryPeakResults results, Rectangle bounds, int[] limits) {
    ArrayList<Spot> list = new ArrayList<Spot>(limits[1] - limits[0] + 1);
    final float minx = bounds.x;
    final float miny = bounds.y;
    final float maxx = bounds.x + bounds.width;
    final float maxy = bounds.y + bounds.height;
    // Find spots within the ROI
    for (PeakResult r : results) {
        final float x = r.params[Gaussian2DFunction.X_POSITION];
        if (x > minx && x < maxx) {
            final float y = r.params[Gaussian2DFunction.Y_POSITION];
            if (y > miny && y < maxy)
                list.add(new Spot(r.getFrame(), x, y, r.getSignal()));
        }
    }
    // For each frame pick the strongest spot
    Collections.sort(list);
    ArrayList<Spot> newList = new ArrayList<Spot>(list.size());
    int currentT = -1;
    for (Spot spot : list) {
        if (currentT != spot.t) {
            newList.add(spot);
            currentT = spot.t;
        }
    }
    return newList.toArray(new Spot[newList.size()]);
}

Example 80 with PeakResult

use of gdsc.smlm.results.PeakResult in project GDSC-SMLM by aherbert.

the class CropResults method cropResults.

/**
	 * Apply the filters to the data
	 */
private void cropResults() {
    MemoryPeakResults newResults = new MemoryPeakResults();
    // These bounds are integer. But this is because the results are meant to come from an image.
    Rectangle bounds = results.getBounds(true);
    // The crop bounds can be floating point...
    // Border
    double xx = bounds.x + border;
    double yy = bounds.y + border;
    double w = Math.max(0, bounds.width - 2 * border);
    double h = Math.max(0, bounds.height - 2 * border);
    Rectangle2D borderBounds = new Rectangle2D.Double(xx, yy, w, h);
    // Bounding box
    if (selectRegion) {
        Rectangle2D boxBounds = new Rectangle2D.Double(x, y, width, height);
        borderBounds = borderBounds.createIntersection(boxBounds);
    }
    // and create another intersection
    if (myUseRoi) {
        ImagePlus imp = WindowManager.getImage(roiImage);
        if (imp != null && imp.getRoi() != null) {
            Rectangle roi = imp.getRoi().getBounds();
            int roiImageWidth = imp.getWidth();
            int roiImageHeight = imp.getHeight();
            double xscale = (double) roiImageWidth / bounds.width;
            double yscale = (double) roiImageHeight / bounds.height;
            Rectangle2D roiBounds = new Rectangle2D.Double(roi.x / xscale, roi.y / yscale, roi.width / xscale, roi.height / yscale);
            borderBounds = borderBounds.createIntersection(roiBounds);
        }
    }
    if (borderBounds.getWidth() > 0 && borderBounds.getHeight() > 0) {
        for (PeakResult result : results.getResults()) {
            if (borderBounds.contains(result.getXPosition(), result.getYPosition()))
                newResults.add(result);
        }
    }
    newResults.copySettings(results);
    newResults.setBounds(new Rectangle((int) Math.floor(borderBounds.getX()), (int) Math.floor(borderBounds.getY()), (int) Math.ceil(borderBounds.getWidth()), (int) Math.ceil(borderBounds.getHeight())));
    if (!overwrite) {
        newResults.setName(results.getName() + " Cropped");
    }
    MemoryPeakResults.addResults(newResults);
    IJ.showStatus(newResults.size() + " Cropped localisations");
}