Search in sources :

Example 11 with FloatProcessor

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

the class FRC method calculateFrcCurve.

/**
	 * Calculate the Fourier Ring Correlation curve for two images.
	 *
	 * @param ip1
	 *            The first image
	 * @param ip2
	 *            The second image
	 * @param nmPerPixel
	 *            the nm per pixel for the super-resolution images
	 * @return the FRC curve
	 */
public FRCCurve calculateFrcCurve(ImageProcessor ip1, ImageProcessor ip2, double nmPerPixel) {
    // Allow a progress tracker to be input
    TrackProgress progess = getTrackProgress();
    progess.incrementProgress(0);
    progess.status("Calculating complex FFT images...");
    // Pad images to the same size if different
    final int maxWidth = FastMath.max(ip1.getWidth(), ip2.getWidth());
    final int maxHeight = FastMath.max(ip1.getHeight(), ip2.getHeight());
    if (FastMath.max(maxWidth, maxHeight) > MAX_SIZE) {
        progess.status("Error calculating FRC curve...");
        progess.incrementProgress(1);
        return null;
    }
    final int fieldOfView = FastMath.max(maxWidth, maxHeight);
    ip1 = pad(ip1, maxWidth, maxHeight);
    ip2 = pad(ip2, maxWidth, maxHeight);
    // The mean of each image after applying the taper
    double mean1, mean2;
    // Real and imaginary components
    float[] re1 = null, im1, re2, im2;
    // Do the first image
    ip1 = getSquareTaperedImage(ip1);
    mean1 = taperedImageMean;
    final int size = ip1.getWidth();
    if (JTRANSFORMS && fourierMethod == FourierMethod.JTRANSFORMS) {
        // Speed up by reusing the FFT object which performs pre-computation
        float[] data = new float[size * size * 2];
        FloatFFT_2D fft = new FloatFFT_2D(size, size);
        // For quadrant swap
        FHT2 fht = new FHT2();
        float[] pixels = (float[]) ip1.getPixels();
        System.arraycopy(pixels, 0, data, 0, pixels.length);
        fft.realForwardFull(data);
        // Get the data
        re1 = pixels;
        im1 = new float[pixels.length];
        for (int i = 0, j = 0; i < data.length; j++) {
            re1[j] = data[i++];
            im1[j] = data[i++];
        }
        fht.swapQuadrants(new FloatProcessor(size, size, re1));
        fht.swapQuadrants(new FloatProcessor(size, size, im1));
        progess.incrementProgress(THIRD);
        ip2 = getSquareTaperedImage(ip2);
        mean2 = taperedImageMean;
        pixels = (float[]) ip2.getPixels();
        System.arraycopy(pixels, 0, data, 0, pixels.length);
        for (int i = pixels.length; i < data.length; i++) data[i] = 0;
        fft.realForwardFull(data);
        // Get the data
        re2 = pixels;
        im2 = new float[pixels.length];
        for (int i = 0, j = 0; i < data.length; j++) {
            re2[j] = data[i++];
            im2[j] = data[i++];
        }
        fht.swapQuadrants(new FloatProcessor(size, size, re2));
        fht.swapQuadrants(new FloatProcessor(size, size, im2));
        progess.incrementProgress(THIRD);
    } else {
        // Simple implementation. This is left for testing.
        //FloatProcessor[] fft = getComplexFFT(ip1);
        //mean1 = taperedImageMean;
        //re1 = (float[]) fft[0].getPixels();
        //im1 = (float[]) fft[1].getPixels();
        //progess.incrementProgress(THIRD);
        //
        //fft = getComplexFFT(ip2);
        //mean2 = taperedImageMean;
        //re2 = (float[]) fft[0].getPixels();
        //im2 = (float[]) fft[1].getPixels();
        //progess.incrementProgress(THIRD);
        // Speed up by reusing the FHT object which performs pre-computation
        FHT2 fht = new FHT2();
        fht.setShowProgress(false);
        float[] f1 = (float[]) ip1.getPixels();
        fht.rc2DFHT(f1, false, size);
        FHT2 fht1 = new FHT2(ip1, true);
        FloatProcessor[] fft = getProcessors(fht1.getComplexTransform2());
        re1 = (float[]) fft[0].getPixels();
        im1 = (float[]) fft[1].getPixels();
        progess.incrementProgress(THIRD);
        ip2 = getSquareTaperedImage(ip2);
        mean2 = taperedImageMean;
        float[] f2 = (float[]) ip2.getPixels();
        fht.rc2DFHT(f2, false, size);
        FHT2 fht2 = new FHT2(ip2, true);
        fft = getProcessors(fht2.getComplexTransform2());
        re2 = (float[]) fft[0].getPixels();
        im2 = (float[]) fft[1].getPixels();
        progess.incrementProgress(THIRD);
    }
    progess.status("Preparing FRC curve calculation...");
    final int centre = size / 2;
    // In-line for speed 
    float[] conjMult = new float[re1.length];
    float[] absFFT1 = new float[re1.length];
    float[] absFFT2 = new float[re1.length];
    // Normalise the FFT to the field of view, i.e. normalise by 1/sqrt(N) for each dimension
    final double norm = 1.0 / fieldOfView;
    for (int i = 0; i < re1.length; i++) {
        re1[i] *= norm;
        im1[i] *= norm;
        re2[i] *= norm;
        im2[i] *= norm;
    }
    boolean basic = false;
    if (basic) {
        compute(conjMult, absFFT1, absFFT2, re1, im1, re2, im2);
    } else {
        computeMirroredFast(size, conjMult, absFFT1, absFFT2, re1, im1, re2, im2);
    }
    progess.status("Calculating FRC curve...");
    final int max = centre - 1;
    FRCCurveResult[] results = new FRCCurveResult[max];
    if (samplingMethod == SamplingMethod.INTERPOLATED_CIRCLE) {
        // Set the results for the centre pixel
        int cx = size * centre + centre;
        results[0] = new FRCCurveResult(0, 1, conjMult[cx], absFFT1[cx], absFFT2[cx]);
        float[][] images = new float[][] { conjMult, absFFT1, absFFT2 };
        for (int radius = 1; radius < max; radius++) {
            // Inline the calculation for speed
            double sum0 = 0;
            double sum1 = 0;
            double sum2 = 0;
            // Note: The image has 2-fold radial symmetry. So we only need to sample
            // angles from 0-pi. To sample the perimeter at pixel intervals we need
            // pi*r samples. So the angle step is max_angle / samples == pi / (pi*r) == 1 / r.
            // The number of samples is increased using the sampling factor.
            final double angleStep = 1 / (perimeterSamplingFactor * radius);
            double angle = 0;
            int numSum = 0;
            while (angle < Math.PI) {
                double cosA = FastMath.cos(angle);
                double x = centre + radius * cosA;
                //double sinA = FastMath.sin(angle);
                double sinA = getSine(angle, cosA);
                double y = centre + radius * sinA;
                double[] values = getInterpolatedValues(x, y, images, size);
                sum0 += values[0];
                sum1 += values[1];
                sum2 += values[2];
                numSum++;
                angle += angleStep;
            }
            results[radius] = new FRCCurveResult(radius, numSum, sum0, sum1, sum2);
        }
    } else {
        // Compute the radial sum as per the DIP image Matlab toolbox
        double[][] sum = RadialStatistics.radialSumAndCount(size, conjMult, absFFT1, absFFT2);
        for (int radius = 0; radius < max; radius++) {
            results[radius] = new FRCCurveResult(radius, (int) sum[3][radius], sum[0][radius], sum[1][radius], sum[2][radius]);
        }
    }
    progess.incrementProgress(LAST_THIRD);
    progess.status("Finished calculating FRC curve...");
    return new FRCCurve(nmPerPixel, fieldOfView, mean1, mean2, results);
}
Also used : FHT2(ij.process.FHT2) FloatProcessor(ij.process.FloatProcessor) NullTrackProgress(gdsc.core.logging.NullTrackProgress) TrackProgress(gdsc.core.logging.TrackProgress) FloatFFT_2D(org.jtransforms.fft.FloatFFT_2D)

Example 12 with FloatProcessor

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

the class IJImagePeakResultsTest method noInterpolateDownInYAtImageEdge.

@Test
public void noInterpolateDownInYAtImageEdge() {
    IJImagePeakResults r = new IJImagePeakResults(title, bounds, 1);
    r.setDisplayFlags(IJImagePeakResults.DISPLAY_WEIGHTED);
    FloatProcessor fp = new FloatProcessor(bounds.width, bounds.height);
    begin(r);
    addValue(r, 1.5f, 0.5f, 2);
    fp.putPixelValue(1, 0, 2);
    r.end();
    float[] expecteds = getImage(fp);
    float[] actuals = getImage(r);
    Assert.assertArrayEquals(expecteds, actuals, 0);
}
Also used : FloatProcessor(ij.process.FloatProcessor) Test(org.junit.Test)

Example 13 with FloatProcessor

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

the class IJImagePeakResultsTest method canInterpolateDownInX.

@Test
public void canInterpolateDownInX() {
    IJImagePeakResults r = new IJImagePeakResults(title, bounds, 1);
    r.setDisplayFlags(IJImagePeakResults.DISPLAY_WEIGHTED);
    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();
    float[] expecteds = getImage(fp);
    float[] actuals = getImage(r);
    Assert.assertArrayEquals(expecteds, actuals, 0);
}
Also used : FloatProcessor(ij.process.FloatProcessor) Test(org.junit.Test)

Example 14 with FloatProcessor

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

the class IJImagePeakResultsTest method canInterpolateDownInYAtPixelEdge.

@Test
public void canInterpolateDownInYAtPixelEdge() {
    IJImagePeakResults r = new IJImagePeakResults(title, bounds, 1);
    r.setDisplayFlags(IJImagePeakResults.DISPLAY_WEIGHTED);
    FloatProcessor fp = new FloatProcessor(bounds.width, bounds.height);
    begin(r);
    addValue(r, 1.5f, 1f, 2);
    fp.putPixelValue(1, 0, 1);
    fp.putPixelValue(1, 1, 1);
    r.end();
    float[] expecteds = getImage(fp);
    float[] actuals = getImage(r);
    Assert.assertArrayEquals(expecteds, actuals, 0);
}
Also used : FloatProcessor(ij.process.FloatProcessor) Test(org.junit.Test)

Example 15 with FloatProcessor

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

the class BinaryDisplay method run.

/* (non-Javadoc)
	 * @see ij.plugin.filter.PlugInFilter#run(ij.process.ImageProcessor)
	 */
public void run(ImageProcessor ip) {
    //		float min = Float.POSITIVE_INFINITY;
    //		for (int i=0; i<ip.getPixelCount(); i++)
    //		{
    //			final float value = ip.getf(i);
    //			if (value == 0)
    //				continue;
    //			if (value < min)
    //				min = value; 
    //		}
    //		ip.setMinAndMax(0, min);
    //		imp.updateAndDraw();
    FloatProcessor fp = new FloatProcessor(ip.getWidth(), ip.getHeight());
    float[] data = (float[]) fp.getPixels();
    for (int i = 0; i < ip.getPixelCount(); i++) {
        final float value = ip.getf(i);
        if (value == 0)
            continue;
        data[i] = 1;
    }
    ip.snapshot();
    ip.setPixels(0, fp);
    ip.setMinAndMax(0, 1);
    imp.updateAndDraw();
}
Also used : FloatProcessor(ij.process.FloatProcessor)

Aggregations

FloatProcessor (ij.process.FloatProcessor)49 Test (org.junit.Test)22 ImageProcessor (ij.process.ImageProcessor)11 Rectangle (java.awt.Rectangle)8 FHT2 (ij.process.FHT2)5 ImageStack (ij.ImageStack)4 Point (java.awt.Point)4 LinkedList (java.util.LinkedList)4 Future (java.util.concurrent.Future)4 ClusterPoint (gdsc.core.clustering.ClusterPoint)3 IJImagePeakResults (gdsc.smlm.ij.results.IJImagePeakResults)3 ImagePlus (ij.ImagePlus)3 RandomGenerator (org.apache.commons.math3.random.RandomGenerator)3 Well19937c (org.apache.commons.math3.random.Well19937c)3 BasePoint (gdsc.core.match.BasePoint)2 MaximaSpotFilter (gdsc.smlm.filters.MaximaSpotFilter)2 IJTablePeakResults (gdsc.smlm.ij.results.IJTablePeakResults)2 PeakResult (gdsc.smlm.results.PeakResult)2 ColorProcessor (ij.process.ColorProcessor)2 ShortProcessor (ij.process.ShortProcessor)2