use of org.apache.commons.math3.analysis.integration.SimpsonIntegrator in project GDSC-SMLM by aherbert.
the class JumpDistanceAnalysisTest method canIntegrateProbabilityToCumulativeWithSinglePopulation.
// Commented out as this test always passes
//@Test
public void canIntegrateProbabilityToCumulativeWithSinglePopulation() {
JumpDistanceAnalysis jd = new JumpDistanceAnalysis();
jd.setMinD(0);
jd.setMinFraction(0);
SimpsonIntegrator si = new SimpsonIntegrator(1e-3, 1e-8, 2, SimpsonIntegrator.SIMPSON_MAX_ITERATIONS_COUNT);
for (double d : D) {
final double[] params = new double[] { d };
final JumpDistanceFunction fp = jd.new JumpDistanceFunction(null, d);
JumpDistanceCumulFunction fc = jd.new JumpDistanceCumulFunction(null, null, d);
double x = d / 8;
UnivariateFunction func = new UnivariateFunction() {
public double value(double x) {
return fp.evaluate(x, params);
}
};
for (int i = 1; i < 10; i++, x *= 2) {
double e = fc.evaluate(x, params);
// Integrate
double o = si.integrate(10000, func, 0, x);
//log("Integrate d=%.1f : x=%.1f, e=%f, o=%f, iter=%d, eval=%d\n", d, x, e, o, si.getIterations(),
// si.getEvaluations());
Assert.assertEquals("Failed to integrate", e, o, e * 1e-2);
}
}
}
use of org.apache.commons.math3.analysis.integration.SimpsonIntegrator in project gatk by broadinstitute.
the class CNLOHCaller method calcEffectiveAlpha.
protected double calcEffectiveAlpha(final double[] effectivePhis) {
final QAlphaUnivariateFunction qAlpha = new QAlphaUnivariateFunction(effectivePhis);
final AQAlphaUnivariateFunction aQAlpha = new AQAlphaUnivariateFunction(effectivePhis);
final BaseAbstractUnivariateIntegrator integrator = new SimpsonIntegrator();
final double numerator = integrator.integrate(10, aQAlpha, MIN_E_ELPHA_INTEGRATION_RANGE, MAX_E_ELPHA_INTEGRATION_RANGE);
final double denominator = integrator.integrate(10, qAlpha, MIN_E_ELPHA_INTEGRATION_RANGE, MAX_E_ELPHA_INTEGRATION_RANGE);
return numerator / denominator;
}
use of org.apache.commons.math3.analysis.integration.SimpsonIntegrator in project GDSC-SMLM by aherbert.
the class AiryPSFModel method createAiryDistribution.
private static synchronized void createAiryDistribution() {
if (spline != null)
return;
final double relativeAccuracy = 1e-4;
final double absoluteAccuracy = 1e-8;
final int minimalIterationCount = 3;
final int maximalIterationCount = 32;
UnivariateIntegrator integrator = new SimpsonIntegrator(relativeAccuracy, absoluteAccuracy, minimalIterationCount, maximalIterationCount);
UnivariateFunction f = new UnivariateFunction() {
public double value(double x) {
//return AiryPattern.intensity(x) * 2 * Math.PI * x / (4 * Math.PI);
return AiryPattern.intensity(x) * 0.5 * x;
}
};
// Integrate up to a set number of dark rings
int samples = 1000;
final double step = RINGS[SAMPLE_RINGS] / samples;
double to = 0, from = 0;
r = new double[samples + 1];
sum = new double[samples + 1];
for (int i = 1; i < sum.length; i++) {
from = to;
r[i] = to = step * i;
sum[i] = integrator.integrate(2000, f, from, to) + sum[i - 1];
}
if (DoubleEquality.relativeError(sum[samples], POWER[SAMPLE_RINGS]) > 1e-3)
throw new RuntimeException("Failed to create the Airy distribution");
SplineInterpolator si = new SplineInterpolator();
spline = si.interpolate(sum, r);
}
use of org.apache.commons.math3.analysis.integration.SimpsonIntegrator in project GDSC-SMLM by aherbert.
the class PoissonFunctionTest method cumulativeProbabilityIsOneWithRealAbove4.
private void cumulativeProbabilityIsOneWithRealAbove4(final double gain, final double mu, int min, int max) {
final double o = mu * gain;
final PoissonFunction f = new PoissonFunction(1.0 / gain, true);
double p = 0;
SimpsonIntegrator in = new SimpsonIntegrator(3, 30);
try {
p = in.integrate(Integer.MAX_VALUE, new UnivariateFunction() {
public double value(double x) {
return f.likelihood(x, o);
}
}, min, max);
System.out.printf("g=%f, mu=%f, o=%f, p=%f\n", gain, mu, o, p);
Assert.assertEquals(String.format("g=%f, mu=%f", gain, mu), 1, p, 0.02);
} catch (TooManyEvaluationsException e) {
//double inc = max / 20000.0;
//for (double x = 0; x <= max; x += inc)
//{
// final double pp = f.likelihood(x, o);
// //System.out.printf("g=%f, mu=%f, o=%f, p=%f\n", gain, mu, o, pp);
// p += pp;
//}
//System.out.printf("g=%f, mu=%f, o=%f, p=%f\n", gain, mu, o, p);
Assert.assertFalse(e.getMessage(), true);
}
}
use of org.apache.commons.math3.analysis.integration.SimpsonIntegrator in project GDSC-SMLM by aherbert.
the class PoissonCalculatorTest method cumulativeProbabilityIsOneWithRealData.
private void cumulativeProbabilityIsOneWithRealData(final double mu, double min, double max, boolean test) {
double p = 0;
UnivariateIntegrator in = new SimpsonIntegrator();
p = in.integrate(20000, new UnivariateFunction() {
public double value(double x) {
double v;
v = PoissonCalculator.likelihood(mu, x);
//System.out.printf("x=%f, v=%f\n", x, v);
return v;
}
}, min, max);
System.out.printf("mu=%f, p=%f\n", mu, p);
if (test) {
Assert.assertEquals(String.format("mu=%f", mu), P_LIMIT, p, 0.02);
}
}
Aggregations