Search in sources :

Example 6 with IndexSupplier

use of uk.ac.sussex.gdsc.test.utils.functions.IndexSupplier in project GDSC-SMLM by aherbert.

the class LsqLvmGradientProcedureTest method gradientProcedureLinearIsFasterThanGradientProcedureMatrix.

private void gradientProcedureLinearIsFasterThanGradientProcedureMatrix(RandomSeed seed, final int nparams, final BaseLsqLvmGradientProcedureFactory factory1, final BaseLsqLvmGradientProcedureFactory factory2, boolean doAssert) {
    Assumptions.assumeTrue(TestSettings.allow(TestComplexity.MEDIUM));
    final int iter = 100;
    final ArrayList<double[]> paramsList = new ArrayList<>(iter);
    final ArrayList<double[]> yList = new ArrayList<>(iter);
    createData(RngUtils.create(seed.getSeed()), 1, iter, paramsList, yList);
    // Remove the timing of the function call by creating a dummy function
    final Gradient1Function func = new FakeGradientFunction(blockWidth, nparams);
    // Create messages
    final IndexSupplier msgA = new IndexSupplier(1, "A ", null);
    final IndexSupplier msgB = new IndexSupplier(1, "B ", null);
    for (int i = 0; i < paramsList.size(); i++) {
        final BaseLsqLvmGradientProcedure p1 = factory1.createProcedure(yList.get(i), func);
        p1.gradient(paramsList.get(i));
        p1.gradient(paramsList.get(i));
        final BaseLsqLvmGradientProcedure p2 = factory2.createProcedure(yList.get(i), func);
        p2.gradient(paramsList.get(i));
        p2.gradient(paramsList.get(i));
        // Check they are the same
        Assertions.assertArrayEquals(p1.getAlphaLinear(), p2.getAlphaLinear(), msgA.set(0, i));
        Assertions.assertArrayEquals(p1.beta, p2.beta, msgB.set(0, i));
    }
    // Realistic loops for an optimisation
    final int loops = 15;
    // Run till stable timing
    final Timer t1 = new Timer() {

        @Override
        void run() {
            for (int i = 0, k = 0; i < paramsList.size(); i++) {
                final BaseLsqLvmGradientProcedure p = factory1.createProcedure(yList.get(i), func);
                for (int j = loops; j-- > 0; ) {
                    p.gradient(paramsList.get(k++ % iter));
                }
            }
        }
    };
    final long time1 = t1.getTime();
    final Timer t2 = new Timer(t1.loops) {

        @Override
        void run() {
            for (int i = 0, k = 0; i < paramsList.size(); i++) {
                final BaseLsqLvmGradientProcedure p2 = factory2.createProcedure(yList.get(i), func);
                for (int j = loops; j-- > 0; ) {
                    p2.gradient(paramsList.get(k++ % iter));
                }
            }
        }
    };
    final long time2 = t2.getTime();
    final LogRecord record = TestLogUtils.getTimingRecord(factory1.getClass().getSimpleName() + nparams, time1, factory2.getClass().getSimpleName(), time2);
    if (!doAssert && record.getLevel() == TestLevel.TEST_FAILURE) {
        record.setLevel(TestLevel.TEST_WARNING);
    }
    logger.log(record);
}
Also used : Gradient1Function(uk.ac.sussex.gdsc.smlm.function.Gradient1Function) IndexSupplier(uk.ac.sussex.gdsc.test.utils.functions.IndexSupplier) LogRecord(java.util.logging.LogRecord) ArrayList(java.util.ArrayList) FakeGradientFunction(uk.ac.sussex.gdsc.smlm.function.FakeGradientFunction)

Example 7 with IndexSupplier

use of uk.ac.sussex.gdsc.test.utils.functions.IndexSupplier in project GDSC-SMLM by aherbert.

the class LvmGradientProcedureTest method gradientProcedureUnrolledComputesSameAsGradientProcedure.

private void gradientProcedureUnrolledComputesSameAsGradientProcedure(RandomSeed seed, int nparams, Type type, boolean precomputed) {
    final int iter = 10;
    final ArrayList<double[]> paramsList = new ArrayList<>(iter);
    final ArrayList<double[]> yList = new ArrayList<>(iter);
    createFakeData(RngUtils.create(seed.getSeed()), nparams, iter, paramsList, yList);
    Gradient1Function func = new FakeGradientFunction(blockWidth, nparams);
    if (precomputed) {
        final double[] b = SimpleArrayUtils.newArray(func.size(), 0.1, 1.3);
        func = OffsetGradient1Function.wrapGradient1Function(func, b);
    }
    final FastLog fastLog = type == Type.FAST_LOG_MLE ? getFastLog() : null;
    final String name = String.format("[%d] %b", nparams, type);
    // Create messages
    final IndexSupplier msgR = new IndexSupplier(1, name + "Result: Not same ", null);
    final IndexSupplier msgOb = new IndexSupplier(1, name + "Observations: Not same beta ", null);
    final IndexSupplier msgOal = new IndexSupplier(1, name + "Observations: Not same alpha linear ", null);
    final IndexSupplier msgOam = new IndexSupplier(1, name + "Observations: Not same alpha matrix ", null);
    for (int i = 0; i < paramsList.size(); i++) {
        final LvmGradientProcedure p1 = createProcedure(type, yList.get(i), func, fastLog);
        p1.gradient(paramsList.get(i));
        final LvmGradientProcedure p2 = LvmGradientProcedureUtils.create(yList.get(i), func, type, fastLog);
        p2.gradient(paramsList.get(i));
        // Exactly the same ...
        Assertions.assertEquals(p1.value, p2.value, msgR.set(0, i));
        Assertions.assertArrayEquals(p1.beta, p2.beta, msgOb.set(0, i));
        Assertions.assertArrayEquals(p1.getAlphaLinear(), p2.getAlphaLinear(), msgOal.set(0, i));
        final double[][] am1 = p1.getAlphaMatrix();
        final double[][] am2 = p2.getAlphaMatrix();
        Assertions.assertArrayEquals(am1, am2, msgOam.set(0, i));
    }
}
Also used : Gradient1Function(uk.ac.sussex.gdsc.smlm.function.Gradient1Function) OffsetGradient1Function(uk.ac.sussex.gdsc.smlm.function.OffsetGradient1Function) IndexSupplier(uk.ac.sussex.gdsc.test.utils.functions.IndexSupplier) ArrayList(java.util.ArrayList) FakeGradientFunction(uk.ac.sussex.gdsc.smlm.function.FakeGradientFunction) FastLog(uk.ac.sussex.gdsc.smlm.function.FastLog)

Example 8 with IndexSupplier

use of uk.ac.sussex.gdsc.test.utils.functions.IndexSupplier in project GDSC-SMLM by aherbert.

the class ParameterBoundsTest method canStepParameter.

private static void canStepParameter(double value) {
    final ParameterBounds bounds = new ParameterBounds(new FakeGradientFunction(1, 1, 1));
    double[] a1 = new double[1];
    double[] a2 = new double[1];
    double[] tmp;
    final double[] step = new double[] { value };
    final IndexSupplier msg = new IndexSupplier(1, "Step ", null);
    for (int i = 1; i <= 10; i++) {
        bounds.applyBounds(a1, step, a2);
        Assertions.assertArrayEquals(a2, new double[] { i * value }, msg.set(0, i));
        tmp = a1;
        a1 = a2;
        a2 = tmp;
    }
}
Also used : IndexSupplier(uk.ac.sussex.gdsc.test.utils.functions.IndexSupplier) FakeGradientFunction(uk.ac.sussex.gdsc.smlm.function.FakeGradientFunction)

Example 9 with IndexSupplier

use of uk.ac.sussex.gdsc.test.utils.functions.IndexSupplier in project GDSC-SMLM by aherbert.

the class FhtFilterTest method canFilter.

private static void canFilter(RandomSeed seed, Operation operation) {
    final int size = 16;
    final int ex = 5;
    final int ey = 7;
    final int ox = 1;
    final int oy = 2;
    final UniformRandomProvider r = RngUtils.create(seed.getSeed());
    final FloatProcessor fp1 = createProcessor(size, ex, ey, 4, 4, r);
    // This is offset from the centre
    final FloatProcessor fp2 = createProcessor(size, size / 2 + ox, size / 2 + oy, 4, 4, r);
    final float[] input1 = ((float[]) fp1.getPixels()).clone();
    final float[] input2 = ((float[]) fp2.getPixels()).clone();
    final FHT fht1 = new FHT(fp1);
    fht1.transform();
    final FHT fht2 = new FHT(fp2);
    fht2.transform();
    FHT fhtE;
    switch(operation) {
        case CONVOLUTION:
            fhtE = fht1.multiply(fht2);
            break;
        case CORRELATION:
            fhtE = fht1.conjugateMultiply(fht2);
            break;
        case DECONVOLUTION:
            fhtE = fht1.divide(fht2);
            break;
        default:
            throw new RuntimeException();
    }
    fhtE.inverseTransform();
    fhtE.swapQuadrants();
    final float[] e = (float[]) fhtE.getPixels();
    if (operation == Operation.CORRELATION) {
        // Test the max correlation position
        final int max = SimpleArrayUtils.findMaxIndex(e);
        final int x = max % 16;
        final int y = max / 16;
        Assertions.assertEquals(ex, x + ox);
        Assertions.assertEquals(ey, y + oy);
    }
    // Test verses a spatial domain filter in the middle of the image
    if (operation != Operation.DECONVOLUTION) {
        double sum = 0;
        float[] i2 = input2;
        if (operation == Operation.CONVOLUTION) {
            i2 = i2.clone();
            KernelFilter.rotate180(i2);
        }
        for (int i = 0; i < input1.length; i++) {
            sum += input1[i] * i2[i];
        }
        // double exp = e[size / 2 * size + size / 2];
        // logger.fine(() -> String.format("Sum = %f vs [%d] %f", sum, size / 2 * size + size / 2,
        // exp);
        Assertions.assertEquals(sum, sum, 1e-3);
    }
    // Test the FHT filter
    final FhtFilter ff = new FhtFilter(input2, size, size);
    ff.setOperation(operation);
    ff.filter(input1, size, size);
    // There may be differences due to the use of the JTransforms library
    final double error = (operation == Operation.DECONVOLUTION) ? 5e-2 : 1e-4;
    final FloatFloatBiPredicate predicate = TestHelper.floatsAreClose(error, 0);
    // This tests everything and can fail easily depending on the random generator
    // due to edge artifacts.
    // TestAssertions.assertArrayTest(e, input1, TestHelper.almostEqualFloats(error, 0));
    // This tests the centre to ignore edge differences
    final int min = size / 4;
    final int max = size - min;
    int repeats = 0;
    for (int y = min; y < max; y++) {
        for (int x = min; x < max; x++) {
            repeats++;
        }
    }
    // Use a fail counter for a 'soft' test that detects major problems
    final int failureLimit = TestCounter.computeFailureLimit(repeats, 0.1);
    final TestCounter failCounter = new TestCounter(failureLimit);
    final IndexSupplier msg = new IndexSupplier(2);
    for (int y = min; y < max; y++) {
        msg.set(1, y);
        for (int x = min; x < max; x++) {
            final int xx = x;
            final int i = y * size + x;
            failCounter.run(() -> {
                TestAssertions.assertTest(e[i], input1[i], predicate, msg.set(0, xx));
            });
        }
    }
}
Also used : TestCounter(uk.ac.sussex.gdsc.test.utils.TestCounter) FHT(ij.process.FHT) FloatProcessor(ij.process.FloatProcessor) IndexSupplier(uk.ac.sussex.gdsc.test.utils.functions.IndexSupplier) FloatFloatBiPredicate(uk.ac.sussex.gdsc.test.api.function.FloatFloatBiPredicate) UniformRandomProvider(org.apache.commons.rng.UniformRandomProvider) FhtFilter(uk.ac.sussex.gdsc.smlm.ij.filters.FhtFilter)

Aggregations

IndexSupplier (uk.ac.sussex.gdsc.test.utils.functions.IndexSupplier)9 FakeGradientFunction (uk.ac.sussex.gdsc.smlm.function.FakeGradientFunction)8 ArrayList (java.util.ArrayList)7 Gradient1Function (uk.ac.sussex.gdsc.smlm.function.Gradient1Function)3 DenseMatrix64F (org.ejml.data.DenseMatrix64F)2 FastLog (uk.ac.sussex.gdsc.smlm.function.FastLog)2 OffsetGradient1Function (uk.ac.sussex.gdsc.smlm.function.OffsetGradient1Function)2 FHT (ij.process.FHT)1 FloatProcessor (ij.process.FloatProcessor)1 LogRecord (java.util.logging.LogRecord)1 UniformRandomProvider (org.apache.commons.rng.UniformRandomProvider)1 FhtFilter (uk.ac.sussex.gdsc.smlm.ij.filters.FhtFilter)1 DoubleDoubleBiPredicate (uk.ac.sussex.gdsc.test.api.function.DoubleDoubleBiPredicate)1 FloatFloatBiPredicate (uk.ac.sussex.gdsc.test.api.function.FloatFloatBiPredicate)1 TestCounter (uk.ac.sussex.gdsc.test.utils.TestCounter)1