Search in sources :

Example 16 with Gaussian2DFunction

use of uk.ac.sussex.gdsc.smlm.function.gaussian.Gaussian2DFunction in project GDSC-SMLM by aherbert.

the class ErfGaussian2DFunctionTest method factoryDefaultsToErfGaussian2DFunction.

@Test
void factoryDefaultsToErfGaussian2DFunction() {
    Gaussian2DFunction f1;
    Gaussian2DFunction f2;
    final int flags2 = BitFlagUtils.unset(flags, GaussianFunctionFactory.FIT_ERF);
    if (this.f2 != null) {
        f1 = GaussianFunctionFactory.create2D(2, maxx, maxy, flags, zModel);
        f2 = GaussianFunctionFactory.create2D(2, maxx, maxy, flags2, zModel);
        Assertions.assertTrue(f1.getClass() == f2.getClass(), "Incorrect function2");
    } else {
        f1 = GaussianFunctionFactory.create2D(1, maxx, maxy, flags, zModel);
        f2 = GaussianFunctionFactory.create2D(1, maxx, maxy, flags2, zModel);
        Assertions.assertTrue(f1.getClass() == f2.getClass(), "Incorrect function1");
    }
}
Also used : Gaussian2DFunction(uk.ac.sussex.gdsc.smlm.function.gaussian.Gaussian2DFunction) Gaussian2DFunctionTest(uk.ac.sussex.gdsc.smlm.function.gaussian.Gaussian2DFunctionTest) Test(org.junit.jupiter.api.Test)

Example 17 with Gaussian2DFunction

use of uk.ac.sussex.gdsc.smlm.function.gaussian.Gaussian2DFunction in project GDSC-SMLM by aherbert.

the class ErfGaussian2DFunctionVsPsfModelTest method computesSameAsPsfModel.

private void computesSameAsPsfModel(double sum, double x0, double x1, double s0, double s1) {
    final Gaussian2DFunction f = GaussianFunctionFactory.create2D(1, width, height, GaussianFunctionFactory.FIT_ERF_FREE_CIRCLE, null);
    final double[] a = new double[1 + Gaussian2DFunction.PARAMETERS_PER_PEAK];
    a[Gaussian2DFunction.SIGNAL] = sum;
    a[Gaussian2DFunction.X_POSITION] = x0;
    a[Gaussian2DFunction.Y_POSITION] = x1;
    a[Gaussian2DFunction.X_SD] = s0;
    a[Gaussian2DFunction.Y_SD] = s1;
    final double[] o = new StandardValueProcedure().getValues(f, a);
    final GaussianPsfModel m = new GaussianPsfModel(s0, s1);
    final double[] e = new double[o.length];
    // Note that the Gaussian2DFunction has 0,0 at the centre of a pixel.
    // The model has 0.5,0.5 at the centre so add an offset.
    m.create2D(e, width, height, sum, x0 + 0.5, x1 + 0.5, null);
    for (int i = 0; i < e.length; i++) {
        // Only check where there is a reasonable amount of signal
        if (e[i] > 1e-2) {
            final double error = DoubleEquality.relativeError(e[i], o[i]);
            // uses the Apache commons implementation.
            if (error > 5e-3) {
                Assertions.fail(String.format("[%d] %s != %s  error = %f", i, Double.toString(e[i]), Double.toString(o[i]), error));
            }
        }
    }
}
Also used : GaussianPsfModel(uk.ac.sussex.gdsc.smlm.model.GaussianPsfModel) Gaussian2DFunction(uk.ac.sussex.gdsc.smlm.function.gaussian.Gaussian2DFunction) StandardValueProcedure(uk.ac.sussex.gdsc.smlm.function.StandardValueProcedure)

Example 18 with Gaussian2DFunction

use of uk.ac.sussex.gdsc.smlm.function.gaussian.Gaussian2DFunction in project GDSC-SMLM by aherbert.

the class PoissonLikelihoodWrapperTest method functionComputesGradient.

private void functionComputesGradient(int flags) {
    final Gaussian2DFunction f1 = GaussianFunctionFactory.create2D(1, maxx, maxx, flags, null);
    testbackground = testbackgroundOptions;
    testsignal1 = testsignal1Options;
    testcx1 = testcx1Options;
    testcy1 = testcy1Options;
    testcz1 = testcz1Options;
    testw1 = testw1Options;
    testangle1 = testangle1Options;
    if (!f1.evaluatesBackground()) {
        testbackground = new double[] { testbackground[0] };
    }
    if (!f1.evaluatesSignal()) {
        testsignal1 = new double[] { testsignal1[0] };
    }
    if (!f1.evaluatesZ()) {
        testcz1 = new double[] { 0 };
    }
    boolean noSecondWidth = false;
    if (!f1.evaluatesSD0()) {
        // Just use 1 width
        testw1 = new double[][] { testw1[0] };
        // If no width 0 then assume we have no width 1 as well
        noSecondWidth = true;
    } else if (!f1.evaluatesSD1()) {
        // No evaluation of second width needs only variation in width 0 so truncate
        testw1 = Arrays.copyOf(testw1, 2);
        noSecondWidth = true;
    }
    if (noSecondWidth) {
        for (int i = 0; i < testw1.length; i++) {
            testw1[i][1] = testw1[i][0];
        }
    }
    if (!f1.evaluatesAngle()) {
        testangle1 = new double[] { 0 };
    }
    final double fraction = 90;
    if (f1.evaluatesBackground()) {
        functionComputesTargetGradient(f1, Gaussian2DFunction.BACKGROUND, fraction);
    }
    if (f1.evaluatesSignal()) {
        functionComputesTargetGradient(f1, Gaussian2DFunction.SIGNAL, fraction);
    }
    functionComputesTargetGradient(f1, Gaussian2DFunction.X_POSITION, fraction);
    functionComputesTargetGradient(f1, Gaussian2DFunction.Y_POSITION, fraction);
    if (f1.evaluatesZ()) {
        functionComputesTargetGradient(f1, Gaussian2DFunction.Z_POSITION, fraction);
    }
    if (f1.evaluatesSD0()) {
        functionComputesTargetGradient(f1, Gaussian2DFunction.X_SD, fraction);
    }
    if (f1.evaluatesSD1()) {
        functionComputesTargetGradient(f1, Gaussian2DFunction.Y_SD, fraction);
    }
    if (f1.evaluatesAngle()) {
        functionComputesTargetGradient(f1, Gaussian2DFunction.ANGLE, fraction);
    }
}
Also used : Gaussian2DFunction(uk.ac.sussex.gdsc.smlm.function.gaussian.Gaussian2DFunction)

Example 19 with Gaussian2DFunction

use of uk.ac.sussex.gdsc.smlm.function.gaussian.Gaussian2DFunction in project GDSC-SMLM by aherbert.

the class ErfGaussian2DFunctionTest method functionIsFasterThanEquivalentGaussian2DFunction.

// Speed test verses equivalent Gaussian2DFunction
@SpeedTag
@Test
void functionIsFasterThanEquivalentGaussian2DFunction() {
    Assumptions.assumeTrue(TestSettings.allow(TestComplexity.MEDIUM));
    final int flags = this.flags & ~GaussianFunctionFactory.FIT_ERF;
    final Gaussian2DFunction gf = GaussianFunctionFactory.create2D(1, maxx, maxy, flags, zModel);
    final boolean zDepth = (flags & GaussianFunctionFactory.FIT_Z) != 0;
    final LocalList<double[]> params1 = new LocalList<>();
    final LocalList<double[]> params2 = new LocalList<>();
    for (final double background : testbackground) {
        // Peak 1
        for (final double signal1 : testsignal1) {
            for (final double cx1 : testcx1) {
                for (final double cy1 : testcy1) {
                    for (final double cz1 : testcz1) {
                        for (final double[] w1 : testw1) {
                            for (final double angle1 : testangle1) {
                                double[] params = createParameters(background, signal1, cx1, cy1, cz1, w1[0], w1[1], angle1);
                                params1.add(params);
                                if (zDepth) {
                                    // Change to a standard free circular function
                                    params = params.clone();
                                    params[Gaussian2DFunction.X_SD] *= zModel.getSx(params[Gaussian2DFunction.Z_POSITION]);
                                    params[Gaussian2DFunction.Y_SD] *= zModel.getSy(params[Gaussian2DFunction.Z_POSITION]);
                                    params[Gaussian2DFunction.Z_POSITION] = 0;
                                    params2.add(params);
                                }
                            }
                        }
                    }
                }
            }
        }
    }
    final double[][] x = params1.toArray(new double[0][]);
    final double[][] x2 = (zDepth) ? params2.toArray(new double[0][]) : x;
    final int runs = 10000 / x.length;
    final TimingService ts = new TimingService(runs);
    ts.execute(new FunctionTimingTask(gf, x2, 1));
    ts.execute(new FunctionTimingTask(gf, x2, 0));
    ts.execute(new FunctionTimingTask(f1, x, 2));
    ts.execute(new FunctionTimingTask(f1, x, 1));
    ts.execute(new FunctionTimingTask(f1, x, 0));
    final int size = ts.getSize();
    ts.repeat(size);
    if (logger.isLoggable(Level.INFO)) {
        logger.info(ts.getReport());
    }
    for (int i = 1; i <= 2; i++) {
        final TimingResult slow = ts.get(-i - 3);
        final TimingResult fast = ts.get(-i);
        logger.log(TestLogUtils.getTimingRecord(slow, fast));
    }
}
Also used : LocalList(uk.ac.sussex.gdsc.core.utils.LocalList) TimingResult(uk.ac.sussex.gdsc.test.utils.TimingResult) Gaussian2DFunction(uk.ac.sussex.gdsc.smlm.function.gaussian.Gaussian2DFunction) TimingService(uk.ac.sussex.gdsc.test.utils.TimingService) SpeedTag(uk.ac.sussex.gdsc.test.junit5.SpeedTag) Gaussian2DFunctionTest(uk.ac.sussex.gdsc.smlm.function.gaussian.Gaussian2DFunctionTest) Test(org.junit.jupiter.api.Test)

Example 20 with Gaussian2DFunction

use of uk.ac.sussex.gdsc.smlm.function.gaussian.Gaussian2DFunction in project GDSC-SMLM by aherbert.

the class ScmosLikelihoodWrapperTest method functionComputesGradientPerDatum.

private void functionComputesGradientPerDatum(RandomSeed seed, int flags) {
    final Gaussian2DFunction f1 = GaussianFunctionFactory.create2D(1, maxx, maxx, flags, null);
    // Setup
    testbackground = testbackgroundOptions;
    testsignal1 = testsignal1Options;
    testcx1 = testcx1Options;
    testcy1 = testcy1Options;
    testcz1 = testcz1Options;
    testw1 = testw1Options;
    testangle1 = testangle1Options;
    if (!f1.evaluatesBackground()) {
        testbackground = new double[] { testbackground[0] };
    }
    if (!f1.evaluatesSignal()) {
        testsignal1 = new double[] { testsignal1[0] };
    }
    if (!f1.evaluatesZ()) {
        testcz1 = new double[] { 0 };
    }
    boolean noSecondWidth = false;
    if (!f1.evaluatesSD0()) {
        // Just use 1 width
        testw1 = new double[][] { testw1[0] };
        // If no width 0 then assume we have no width 1 as well
        noSecondWidth = true;
    } else if (!f1.evaluatesSD1()) {
        // No evaluation of second width needs only variation in width 0 so truncate
        testw1 = Arrays.copyOf(testw1, 2);
        noSecondWidth = true;
    }
    if (noSecondWidth) {
        for (int i = 0; i < testw1.length; i++) {
            testw1[i][1] = testw1[i][0];
        }
    }
    if (!f1.evaluatesAngle()) {
        testangle1 = new double[] { 0 };
    }
    if (f1.evaluatesBackground()) {
        functionComputesTargetGradientPerDatum(seed, f1, Gaussian2DFunction.BACKGROUND);
    }
    if (f1.evaluatesSignal()) {
        functionComputesTargetGradientPerDatum(seed, f1, Gaussian2DFunction.SIGNAL);
    }
    functionComputesTargetGradientPerDatum(seed, f1, Gaussian2DFunction.X_POSITION);
    functionComputesTargetGradientPerDatum(seed, f1, Gaussian2DFunction.Y_POSITION);
    if (f1.evaluatesZ()) {
        functionComputesTargetGradientPerDatum(seed, f1, Gaussian2DFunction.Z_POSITION);
    }
    if (f1.evaluatesSD0()) {
        functionComputesTargetGradientPerDatum(seed, f1, Gaussian2DFunction.X_SD);
    }
    if (f1.evaluatesSD1()) {
        functionComputesTargetGradientPerDatum(seed, f1, Gaussian2DFunction.Y_SD);
    }
    if (f1.evaluatesAngle()) {
        functionComputesTargetGradientPerDatum(seed, f1, Gaussian2DFunction.ANGLE);
    }
}
Also used : Gaussian2DFunction(uk.ac.sussex.gdsc.smlm.function.gaussian.Gaussian2DFunction)

Aggregations

Gaussian2DFunction (uk.ac.sussex.gdsc.smlm.function.gaussian.Gaussian2DFunction)33 ErfGaussian2DFunction (uk.ac.sussex.gdsc.smlm.function.gaussian.erf.ErfGaussian2DFunction)7 UniformRandomProvider (org.apache.commons.rng.UniformRandomProvider)6 StandardValueProcedure (uk.ac.sussex.gdsc.smlm.function.StandardValueProcedure)5 SeededTest (uk.ac.sussex.gdsc.test.junit5.SeededTest)5 Test (org.junit.jupiter.api.Test)4 LocalList (uk.ac.sussex.gdsc.core.utils.LocalList)4 StandardFloatValueProcedure (uk.ac.sussex.gdsc.smlm.function.StandardFloatValueProcedure)4 TimingService (uk.ac.sussex.gdsc.test.utils.TimingService)4 ArrayList (java.util.ArrayList)3 GradientCalculator (uk.ac.sussex.gdsc.smlm.fitting.nonlinear.gradient.GradientCalculator)3 ValueProcedure (uk.ac.sussex.gdsc.smlm.function.ValueProcedure)3 DoubleDoubleBiPredicate (uk.ac.sussex.gdsc.test.api.function.DoubleDoubleBiPredicate)3 DenseMatrix64F (org.ejml.data.DenseMatrix64F)2 DoubleEquality (uk.ac.sussex.gdsc.core.utils.DoubleEquality)2 Statistics (uk.ac.sussex.gdsc.core.utils.Statistics)2 PoissonGradientProcedure (uk.ac.sussex.gdsc.smlm.fitting.nonlinear.gradient.PoissonGradientProcedure)2 FisherInformation (uk.ac.sussex.gdsc.smlm.function.FisherInformation)2 Gradient1Function (uk.ac.sussex.gdsc.smlm.function.Gradient1Function)2 Gradient2Function (uk.ac.sussex.gdsc.smlm.function.Gradient2Function)2