Examples with FloatProcessor ij.process.FloatProcessor used on opensource projects

Example 71 with FloatProcessor

use of ij.process.FloatProcessor in project GDSC-SMLM by aherbert.

the class ImageJImagePeakResultsTest method canInterpolateDownInX.

@Test
void canInterpolateDownInX() {
    final ImageJImagePeakResults r = new ImageJImagePeakResults(title, bounds, 1);
    r.setDisplayFlags(ImageJImagePeakResults.DISPLAY_WEIGHTED);
    final FloatProcessor fp = new FloatProcessor(bounds.width, bounds.height);
    begin(r);
    addValue(r, 1.25f, 1.5f, 2);
    fp.putPixelValue(0, 1, 0.5f);
    fp.putPixelValue(1, 1, 1.5f);
    r.end();
    final float[] expecteds = getImage(fp);
    final float[] actuals = getImage(r);
    Assertions.assertArrayEquals(expecteds, actuals);
}

Example 72 with FloatProcessor

use of ij.process.FloatProcessor in project GDSC-SMLM by aherbert.

the class PsfCombiner method combineImages.

private void combineImages() {
    final double nmPerPixel = getNmPerPixel();
    if (nmPerPixel <= 0) {
        return;
    }
    final double nmPerSlice = getNmPerSlice();
    if (nmPerPixel <= 0) {
        return;
    }
    // Find the lowest & highest dimensions
    int minStart = Integer.MAX_VALUE;
    int maxStart = Integer.MIN_VALUE;
    int minEnd = Integer.MAX_VALUE;
    int maxEnd = Integer.MIN_VALUE;
    int minSize = Integer.MAX_VALUE;
    int maxSize = 0;
    for (final Psf psf : input) {
        if (minStart > psf.start) {
            minStart = psf.start;
        }
        if (maxStart < psf.start) {
            maxStart = psf.start;
        }
        if (maxEnd < psf.getEnd()) {
            maxEnd = psf.getEnd();
        }
        if (minEnd > psf.getEnd()) {
            minEnd = psf.getEnd();
        }
        if (maxSize < psf.getSize()) {
            maxSize = psf.getSize();
        }
        if (minSize > psf.getSize()) {
            minSize = psf.getSize();
        }
    }
    int size = maxSize;
    int shift = -minStart;
    int depth = maxEnd - minStart + 1;
    // Option to crop. Do this before processing as it will make the plugin faster
    if (minStart < maxStart || minEnd < maxEnd || minSize < maxSize) {
        boolean crop;
        if (ImageJUtils.isMacro()) {
            final String options = Macro.getOptions();
            crop = options.contains(" crop");
        } else {
            final GenericDialog gd = new GenericDialog(TITLE);
            ImageJUtils.addMessage(gd, "The range of the PSFs is different:\nStart %d to %d\nEnd %d to %d\n" + "Size %d to %d\n \nCrop to the smallest?", minStart, maxStart, minEnd, maxEnd, minSize, maxSize);
            gd.enableYesNoCancel();
            gd.addHelp(HelpUrls.getUrl("psf-combiner"));
            gd.showDialog();
            if (gd.wasCanceled()) {
                return;
            }
            crop = gd.wasOKed();
        }
        if (crop) {
            Recorder.recordOption("crop");
            for (final Psf psf : input) {
                psf.crop(maxStart, minEnd, minSize);
            }
            size = minSize;
            shift = -maxStart;
            depth = minEnd - maxStart + 1;
        }
    }
    // Shift all stacks
    int totalImages = 0;
    for (final Psf psf : input) {
        psf.start += shift;
        totalImages += psf.psfSettings.getImageCount();
    }
    // Create a stack to hold all the images
    // Create a stack to hold the sum of the weights
    final ImageStack stack = new ImageStack(size, size, depth);
    final ImageStack stackW = new ImageStack(size, size, depth);
    for (int n = 1; n <= depth; n++) {
        stack.setPixels(new float[size * size], n);
        stackW.setPixels(new float[size * size], n);
    }
    // Insert all the PSFs
    IJ.showStatus("Creating combined image ...");
    int imageNo = 0;
    final double fraction = 1.0 / input.size();
    for (final Psf psf : input) {
        double progress = imageNo * fraction;
        final ImageStack psfStack = psf.psfStack;
        final int w = psf.getSize();
        final int offsetXy = (size - w) / 2;
        final int offsetZ = psf.start;
        final double weight = (1.0 * psf.psfSettings.getImageCount()) / totalImages;
        final FloatProcessor wp = new FloatProcessor(w, w, SimpleArrayUtils.newFloatArray(psfStack.getWidth() * psfStack.getHeight(), (float) weight));
        final double increment = fraction / psfStack.getSize();
        for (int n = 1; n <= psfStack.getSize(); n++) {
            progress += increment;
            IJ.showProgress(progress);
            // Get the data and adjust using the weight
            final float[] psfData = ImageJImageConverter.getData(psfStack.getProcessor(n));
            for (int i = 0; i < psfData.length; i++) {
                psfData[i] *= weight;
            }
            // Insert into the combined PSF
            final int slice = n + offsetZ;
            ImageProcessor ip = stack.getProcessor(slice);
            ip.copyBits(new FloatProcessor(w, w, psfData), offsetXy, offsetXy, Blitter.ADD);
            // Insert the weights
            ip = stackW.getProcessor(slice);
            ip.copyBits(wp, offsetXy, offsetXy, Blitter.ADD);
        }
        imageNo++;
    }
    // Normalise
    for (int n = 1; n <= depth; n++) {
        stack.getProcessor(n).copyBits(stackW.getProcessor(n), 0, 0, Blitter.DIVIDE);
    }
    ImageJUtils.finished();
    final ImagePlus imp = ImageJUtils.display("Combined PSF", stack);
    imp.setSlice(1 + shift);
    imp.resetDisplayRange();
    imp.updateAndDraw();
    final double fwhm = getFwhm();
    imp.setProperty("Info", ImagePsfHelper.toString(ImagePsfHelper.create(imp.getSlice(), nmPerPixel, nmPerSlice, totalImages, fwhm)));
    ImageJUtils.log("%s : z-centre = %d, nm/Pixel = %s, nm/Slice = %s, %d images, FWHM = %s\n", imp.getTitle(), imp.getSlice(), MathUtils.rounded(nmPerPixel), MathUtils.rounded(nmPerSlice), totalImages, MathUtils.rounded(fwhm));
}

Example 73 with FloatProcessor

use of ij.process.FloatProcessor in project GDSC-SMLM by aherbert.

the class PixelFilter method run.

@Override
public void run(ImageProcessor ip) {
    // Compute rolling sums
    final FloatProcessor fp = ip.toFloat(0, null);
    final float[] data = (float[]) ip.toFloat(0, null).getPixels();
    double[] rollingSum = null;
    double[] rollingSumSq = null;
    if (preview && cachedRollingSum != null) {
        rollingSum = cachedRollingSum;
        rollingSumSq = cachedRollingSumSq;
    }
    if (rollingSum == null || rollingSumSq == null) {
        rollingSum = new double[ip.getPixelCount()];
        rollingSumSq = new double[rollingSum.length];
        calculateRollingSums(fp, rollingSum, rollingSumSq);
    }
    int count = 0;
    final int maxx = ip.getWidth();
    final int maxy = ip.getHeight();
    for (int y = 0, i = 0; y < maxy; y++) {
        for (int x = 0; x < maxx; x++, i++) {
            double sum = 0;
            double sumSquares = 0;
            int minU = x - settings.radius - 1;
            final int maxU = Math.min(x + settings.radius, maxx - 1);
            final int maxV = Math.min(y + settings.radius, maxy - 1);
            // Compute sum from rolling sum using:
            // sum(u,v) =
            // + s(u+N,v+N)
            // - s(u-N-1,v+N)
            // - s(u+N,v-N-1)
            // + s(u-N-1,v-N-1)
            // Note:
            // s(u,v) = 0 when either u,v < 0
            // s(u,v) = s(umax,v) when u>umax
            // s(u,v) = s(u,vmax) when v>vmax
            // s(u,v) = s(umax,vmax) when u>umax,v>vmax
            // Likewise for ss
            // + s(u+N-1,v+N-1)
            int index = maxV * maxx + maxU;
            sum += rollingSum[index];
            sumSquares += rollingSumSq[index];
            if (minU >= 0) {
                // - s(u-1,v+N-1)
                index = maxV * maxx + minU;
                sum -= rollingSum[index];
                sumSquares -= rollingSumSq[index];
            }
            int minV = y - settings.radius - 1;
            if (minV >= 0) {
                // - s(u+N-1,v-1)
                index = minV * maxx + maxU;
                sum -= rollingSum[index];
                sumSquares -= rollingSumSq[index];
                if (minU >= 0) {
                    // + s(u-1,v-1)
                    index = minV * maxx + minU;
                    sum += rollingSum[index];
                    sumSquares += rollingSumSq[index];
                }
            }
            // Reset to bounds to calculate the number of pixels
            if (minU < 0) {
                minU = -1;
            }
            if (minV < 0) {
                minV = -1;
            }
            final int n = (maxU - minU) * (maxV - minV);
            if (n < 2) {
                continue;
            }
            // Get the sum of squared differences
            final double residuals = sumSquares - (sum * sum) / n;
            if (residuals > 0.0) {
                final double stdDev = Math.sqrt(residuals / (n - 1.0));
                final double mean = sum / n;
                if (Math.abs(data[i] - mean) / stdDev > settings.error) {
                    ip.setf(i, (float) mean);
                    count++;
                }
            }
        }
    }
    if (preview) {
        cachedRollingSum = rollingSum;
        cachedRollingSumSq = rollingSumSq;
        label.setText("Replaced " + count);
    } else if (pfr != null && count > 0) {
        IJ.log(String.format("Slice %d : Replaced %d pixels", pfr.getSliceNumber(), count));
    }
}

Example 74 with FloatProcessor

use of ij.process.FloatProcessor in project GDSC-SMLM by aherbert.

the class SpotFinderPreview method run.

private void run(ImageProcessor ip, MaximaSpotFilter filter) {
    if (refreshing) {
        return;
    }
    currentSlice = imp.getCurrentSlice();
    final Rectangle bounds = ip.getRoi();
    // Crop to the ROI
    FloatProcessor fp = ip.crop().toFloat(0, null);
    float[] data = (float[]) fp.getPixels();
    final int width = fp.getWidth();
    final int height = fp.getHeight();
    // Store the mean bias and gain of the region data.
    // This is used to correctly overlay the filtered data on the original image.
    double bias = 0;
    double gain = 1;
    boolean adjust = false;
    // Set weights
    final CameraModel cameraModel = fitConfig.getCameraModel();
    if (!(cameraModel instanceof FakePerPixelCameraModel)) {
        // This should be done on the normalised data
        final float[] w = cameraModel.getNormalisedWeights(bounds);
        filter.setWeights(w, width, height);
        data = data.clone();
        if (data.length < ip.getPixelCount()) {
            adjust = true;
            bias = MathUtils.sum(cameraModel.getBias(bounds)) / data.length;
            gain = MathUtils.sum(cameraModel.getGain(bounds)) / data.length;
        }
        cameraModel.removeBiasAndGain(bounds, data);
    }
    final Spot[] spots = filter.rank(data, width, height);
    data = filter.getPreprocessedData();
    final int size = spots.length;
    if (topNScrollBar != null) {
        topNScrollBar.setMaximum(size);
        selectScrollBar.setMaximum(size);
    }
    fp = new FloatProcessor(width, height, data);
    final FloatProcessor out = new FloatProcessor(ip.getWidth(), ip.getHeight());
    out.copyBits(ip, 0, 0, Blitter.COPY);
    if (adjust) {
        fp.multiply(gain);
        fp.add(bias);
    }
    out.insert(fp, bounds.x, bounds.y);
    final double min = fp.getMin();
    final double max = fp.getMax();
    out.setMinAndMax(min, max);
    final Overlay o = new Overlay();
    o.add(new ImageRoi(0, 0, out));
    if (label != null) {
        // Get results for frame
        final Coordinate[] actual = ResultsMatchCalculator.getCoordinates(actualCoordinates, imp.getCurrentSlice());
        final Coordinate[] predicted = new Coordinate[size];
        for (int i = 0; i < size; i++) {
            predicted[i] = new BasePoint(spots[i].x + bounds.x, spots[i].y + bounds.y);
        }
        // Compute assignments
        final LocalList<FractionalAssignment> fractionalAssignments = new LocalList<>(3 * predicted.length);
        final double matchDistance = settings.distance * fitConfig.getInitialPeakStdDev();
        final RampedScore score = RampedScore.of(matchDistance, matchDistance * settings.lowerDistance / 100, false);
        final double dmin = matchDistance * matchDistance;
        final int nActual = actual.length;
        final int nPredicted = predicted.length;
        for (int j = 0; j < nPredicted; j++) {
            // Centre in the middle of the pixel
            final float x = predicted[j].getX() + 0.5f;
            final float y = predicted[j].getY() + 0.5f;
            // Any spots that match
            for (int i = 0; i < nActual; i++) {
                final double dx = (x - actual[i].getX());
                final double dy = (y - actual[i].getY());
                final double d2 = dx * dx + dy * dy;
                if (d2 <= dmin) {
                    final double d = Math.sqrt(d2);
                    final double s = score.score(d);
                    if (s == 0) {
                        continue;
                    }
                    double distance = 1 - s;
                    if (distance == 0) {
                        // In the case of a match below the distance thresholds
                        // the distance will be 0. To distinguish between candidates all below
                        // the thresholds just take the closest.
                        // We know d2 is below dmin so we subtract the delta.
                        distance -= (dmin - d2);
                    }
                    // Store the match
                    fractionalAssignments.add(new ImmutableFractionalAssignment(i, j, distance, s));
                }
            }
        }
        final FractionalAssignment[] assignments = fractionalAssignments.toArray(new FractionalAssignment[0]);
        // Compute matches
        final RankedScoreCalculator calc = RankedScoreCalculator.create(assignments, nActual - 1, nPredicted - 1);
        final boolean save = settings.showTP || settings.showFP;
        final double[] calcScore = calc.score(nPredicted, settings.multipleMatches, save);
        final ClassificationResult result = RankedScoreCalculator.toClassificationResult(calcScore, nActual);
        // Compute AUC and max jaccard (and plot)
        final double[][] curve = RankedScoreCalculator.getPrecisionRecallCurve(assignments, nActual, nPredicted);
        final double[] precision = curve[0];
        final double[] recall = curve[1];
        final double[] jaccard = curve[2];
        final double auc = AucCalculator.auc(precision, recall);
        // Show scores
        final String scoreLabel = String.format("Slice=%d, AUC=%s, R=%s, Max J=%s", imp.getCurrentSlice(), MathUtils.rounded(auc), MathUtils.rounded(result.getRecall()), MathUtils.rounded(MathUtils.maxDefault(0, jaccard)));
        setLabel(scoreLabel);
        // Plot
        String title = TITLE + " Performance";
        Plot plot = new Plot(title, "Spot Rank", "");
        final double[] rank = SimpleArrayUtils.newArray(precision.length, 0, 1.0);
        plot.setLimits(0, nPredicted, 0, 1.05);
        plot.setColor(Color.blue);
        plot.addPoints(rank, precision, Plot.LINE);
        plot.setColor(Color.red);
        plot.addPoints(rank, recall, Plot.LINE);
        plot.setColor(Color.black);
        plot.addPoints(rank, jaccard, Plot.LINE);
        plot.setColor(Color.black);
        plot.addLabel(0, 0, scoreLabel);
        final WindowOrganiser windowOrganiser = new WindowOrganiser();
        ImageJUtils.display(title, plot, 0, windowOrganiser);
        title = TITLE + " Precision-Recall";
        plot = new Plot(title, "Recall", "Precision");
        plot.setLimits(0, 1, 0, 1.05);
        plot.setColor(Color.red);
        plot.addPoints(recall, precision, Plot.LINE);
        plot.drawLine(recall[recall.length - 1], precision[recall.length - 1], recall[recall.length - 1], 0);
        plot.setColor(Color.black);
        plot.addLabel(0, 0, scoreLabel);
        ImageJUtils.display(title, plot, 0, windowOrganiser);
        windowOrganiser.tile();
        // Create Rois for TP and FP
        if (save) {
            final double[] matchScore = RankedScoreCalculator.getMatchScore(calc.getScoredAssignments(), nPredicted);
            int matches = 0;
            for (int i = 0; i < matchScore.length; i++) {
                if (matchScore[i] != 0) {
                    matches++;
                }
            }
            if (settings.showTP) {
                final float[] x = new float[matches];
                final float[] y = new float[x.length];
                int count = 0;
                for (int i = 0; i < matchScore.length; i++) {
                    if (matchScore[i] != 0) {
                        final BasePoint p = (BasePoint) predicted[i];
                        x[count] = p.getX() + 0.5f;
                        y[count] = p.getY() + 0.5f;
                        count++;
                    }
                }
                addRoi(0, o, x, y, count, Color.green);
            }
            if (settings.showFP) {
                final float[] x = new float[nPredicted - matches];
                final float[] y = new float[x.length];
                int count = 0;
                for (int i = 0; i < matchScore.length; i++) {
                    if (matchScore[i] == 0) {
                        final BasePoint p = (BasePoint) predicted[i];
                        x[count] = p.getX() + 0.5f;
                        y[count] = p.getY() + 0.5f;
                        count++;
                    }
                }
                addRoi(0, o, x, y, count, Color.red);
            }
        }
    } else {
        final WindowOrganiser wo = new WindowOrganiser();
        // Option to show the number of neighbours within a set pixel box radius
        final int[] count = spotFilterHelper.countNeighbours(spots, width, height, settings.neighbourRadius);
        // Show as histogram the totals...
        new HistogramPlotBuilder(TITLE, StoredData.create(count), "Neighbours").setIntegerBins(true).setPlotLabel("Radius = " + settings.neighbourRadius).show(wo);
        // TODO - Draw n=0, n=1 on the image overlay
        final LUT lut = LutHelper.createLut(LutColour.FIRE_LIGHT);
        // These are copied by the ROI
        final float[] x = new float[1];
        final float[] y = new float[1];
        // Plot the intensity
        final double[] intensity = new double[size];
        final double[] rank = SimpleArrayUtils.newArray(size, 1, 1.0);
        final int top = (settings.topN > 0) ? settings.topN : size;
        final int size_1 = size - 1;
        for (int i = 0; i < size; i++) {
            intensity[i] = spots[i].intensity;
            if (i < top) {
                x[0] = spots[i].x + bounds.x + 0.5f;
                y[0] = spots[i].y + bounds.y + 0.5f;
                final Color c = LutHelper.getColour(lut, size_1 - i, size);
                addRoi(0, o, x, y, 1, c, 2, 1);
            }
        }
        final String title = TITLE + " Intensity";
        final Plot plot = new Plot(title, "Rank", "Intensity");
        plot.setColor(Color.blue);
        plot.addPoints(rank, intensity, Plot.LINE);
        if (settings.topN > 0 && settings.topN < size) {
            plot.setColor(Color.magenta);
            plot.drawLine(settings.topN, 0, settings.topN, intensity[settings.topN - 1]);
        }
        if (settings.select > 0 && settings.select < size) {
            plot.setColor(Color.yellow);
            final int index = settings.select - 1;
            final double in = intensity[index];
            plot.drawLine(settings.select, 0, settings.select, in);
            x[0] = spots[index].x + bounds.x + 0.5f;
            y[0] = spots[index].y + bounds.y + 0.5f;
            final Color c = LutHelper.getColour(lut, size_1 - settings.select, size);
            addRoi(0, o, x, y, 1, c, 3, 3);
            plot.setColor(Color.black);
            plot.addLabel(0, 0, "Selected spot intensity = " + MathUtils.rounded(in));
        }
        ImageJUtils.display(title, plot, 0, wo);
        wo.tile();
    }
    imp.setOverlay(o);
}

Example 75 with FloatProcessor

use of ij.process.FloatProcessor in project GDSC-SMLM by aherbert.

the class GaussianFit method addToImage.

/**
 * Adds the function to the image at the specified origin.
 *
 * @param ox the x-origin
 * @param oy the y-origin
 * @param renderedImage the rendered image
 * @param params the function parameters
 * @param npeaks the number of peaks
 * @param width the function width
 * @param height the function height
 */
private void addToImage(int ox, int oy, final FloatProcessor renderedImage, double[] params, int npeaks, int width, int height) {
    final EllipticalGaussian2DFunction f = new EllipticalGaussian2DFunction(npeaks, width, height);
    invertParameters(params);
    final FloatProcessor fp = new FloatProcessor(width, height, f.computeValues(params));
    // Insert into a full size image
    renderedImage.copyBits(fp, ox, oy, Blitter.ADD);
}