Search in sources :

Example 1 with CubicSplinePosition

use of uk.ac.sussex.gdsc.core.math.interpolation.CubicSplinePosition in project GDSC-SMLM by aherbert.

the class CubicSplineCalculatorTest method canComputeCoefficientsForGaussianFunction.

@Test
void canComputeCoefficientsForGaussianFunction() {
    final int x = 4;
    final int y = 4;
    final int z = 4;
    final double xscale = 1;
    final double yscale = 0.5;
    final double zscale = 2.0;
    final double[] xval = SimpleArrayUtils.newArray(x, 0, xscale);
    final double[] yval = SimpleArrayUtils.newArray(y, 0, yscale);
    final double[] zval = SimpleArrayUtils.newArray(z, 0, zscale);
    final double[][][] fval = createData(x, y, z, null);
    final CustomTricubicInterpolatingFunction f1 = new CustomTricubicInterpolator().interpolate(xval, yval, zval, fval);
    final double[] exp = new double[64];
    f1.getSplineNode(1, 1, 1).getCoefficients(exp);
    final CustomTricubicFunction f = CustomTricubicFunctionUtils.create(exp);
    final CubicSplinePosition[] s = new CubicSplinePosition[4];
    for (int i = 0; i < 4; i++) {
        s[i] = new CubicSplinePosition((double) i / 3);
    }
    final double[][][] value = new double[4][4][4];
    final double[] b = new double[64];
    int count = 0;
    for (int k = 0; k < 4; k++) {
        for (int j = 0; j < 4; j++) {
            for (int i = 0; i < 4; i++) {
                value[i][j][k] = f.value(s[i], s[j], s[k]);
                b[count++] = value[i][j][k];
            }
        }
    }
    final CubicSplineCalculator calc = new CubicSplineCalculator();
    double[] obs = calc.compute(value);
    Assertions.assertArrayEquals(exp, obs, 1e-6);
    obs = calc.compute(b);
    Assertions.assertArrayEquals(exp, obs, 1e-6);
}
Also used : CustomTricubicInterpolator(uk.ac.sussex.gdsc.core.math.interpolation.CustomTricubicInterpolator) CustomTricubicInterpolatingFunction(uk.ac.sussex.gdsc.core.math.interpolation.CustomTricubicInterpolatingFunction) CubicSplinePosition(uk.ac.sussex.gdsc.core.math.interpolation.CubicSplinePosition) CustomTricubicFunction(uk.ac.sussex.gdsc.core.math.interpolation.CustomTricubicFunction) Test(org.junit.jupiter.api.Test)

Example 2 with CubicSplinePosition

use of uk.ac.sussex.gdsc.core.math.interpolation.CubicSplinePosition 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 3 with CubicSplinePosition

use of uk.ac.sussex.gdsc.core.math.interpolation.CubicSplinePosition in project GDSC-SMLM by aherbert.

the class CubicSplineData method createCubicSplinePosition.

private static CubicSplinePosition[] createCubicSplinePosition(int n) {
    // Use an extra one to have the final x=1 interpolation point.
    final int n1 = n + 1;
    final double step = 1.0 / n;
    final CubicSplinePosition[] s = new CubicSplinePosition[n1];
    for (int x = 0; x < n; x++) {
        s[x] = new CubicSplinePosition(x * step);
    }
    // Final interpolation point must be exactly 1
    s[n] = new CubicSplinePosition(1);
    return s;
}
Also used : CubicSplinePosition(uk.ac.sussex.gdsc.core.math.interpolation.CubicSplinePosition)

Example 4 with CubicSplinePosition

use of uk.ac.sussex.gdsc.core.math.interpolation.CubicSplinePosition in project GDSC-SMLM by aherbert.

the class CubicSplineCalculatorTest method canComputeCoefficientsForDistanceFunction.

@Test
void canComputeCoefficientsForDistanceFunction() {
    final double[] exp = new double[64];
    int count = 0;
    for (int k = 0; k < 4; k++) {
        for (int j = 0; j < 4; j++) {
            for (int i = 0; i < 4; i++) {
                exp[count++] = Math.sqrt(i * i + j * j + k * k);
            }
        }
    }
    final CustomTricubicFunction f = CustomTricubicFunctionUtils.create(exp);
    final CubicSplinePosition[] s = new CubicSplinePosition[4];
    for (int i = 0; i < 4; i++) {
        s[i] = new CubicSplinePosition((double) i / 3);
    }
    final double[][][] value = new double[4][4][4];
    final double[] b = new double[64];
    count = 0;
    for (int k = 0; k < 4; k++) {
        for (int j = 0; j < 4; j++) {
            for (int i = 0; i < 4; i++) {
                value[i][j][k] = f.value(s[i], s[j], s[k]);
                b[count++] = value[i][j][k];
            }
        }
    }
    final CubicSplineCalculator calc = new CubicSplineCalculator();
    double[] obs = calc.compute(value);
    Assertions.assertArrayEquals(exp, obs, 1e-6);
    obs = calc.compute(b);
    Assertions.assertArrayEquals(exp, obs, 1e-6);
}
Also used : CubicSplinePosition(uk.ac.sussex.gdsc.core.math.interpolation.CubicSplinePosition) CustomTricubicFunction(uk.ac.sussex.gdsc.core.math.interpolation.CustomTricubicFunction) Test(org.junit.jupiter.api.Test)

Aggregations

CubicSplinePosition (uk.ac.sussex.gdsc.core.math.interpolation.CubicSplinePosition)4 CustomTricubicFunction (uk.ac.sussex.gdsc.core.math.interpolation.CustomTricubicFunction)3 Test (org.junit.jupiter.api.Test)2 Ticker (uk.ac.sussex.gdsc.core.logging.Ticker)1 CustomTricubicInterpolatingFunction (uk.ac.sussex.gdsc.core.math.interpolation.CustomTricubicInterpolatingFunction)1 CustomTricubicInterpolator (uk.ac.sussex.gdsc.core.math.interpolation.CustomTricubicInterpolator)1 FloatCustomTricubicFunction (uk.ac.sussex.gdsc.core.math.interpolation.FloatCustomTricubicFunction)1