diff --git a/numass-main/src/main/java/inr/numass/NumassPlugin.java b/numass-main/src/main/java/inr/numass/NumassPlugin.java index 561ff11a..8dff6cc1 100644 --- a/numass-main/src/main/java/inr/numass/NumassPlugin.java +++ b/numass-main/src/main/java/inr/numass/NumassPlugin.java @@ -40,7 +40,6 @@ import org.apache.commons.math3.analysis.UnivariateFunction; import org.apache.commons.math3.util.FastMath; /** - * * @author Alexander Nozik */ @PluginDef(group = "inr.numass", name = "numass", @@ -123,6 +122,18 @@ public class NumassPlugin extends BasicPlugin { double beta = meta.getDouble("tailBeta", 0); return (double E, double U) -> 1 - (E - U) * (alpha + E / 1000d * beta) / 1000d; }); + + math.registerBivariate("numass.resolutionTail.2017", meta -> + (double E, double U) -> { + double D = E - U; + return 0.99797 - 3.05346E-7 * D - 5.45738E-10 * Math.pow(D, 2) - 6.36105E-14 * Math.pow(D, 3); + }); + + math.registerBivariate("numass.resolutionTail.2017.mod", meta -> + (double E, double U) -> { + double D = E - U; + return (0.99797 - 3.05346E-7 * D - 5.45738E-10 * Math.pow(D, 2) - 6.36105E-14 * Math.pow(D, 3)); + }); } /** diff --git a/numass-main/src/main/java/inr/numass/models/LossCalculator.java b/numass-main/src/main/java/inr/numass/models/LossCalculator.java index fe25de38..45ba0ac6 100644 --- a/numass-main/src/main/java/inr/numass/models/LossCalculator.java +++ b/numass-main/src/main/java/inr/numass/models/LossCalculator.java @@ -20,6 +20,7 @@ import hep.dataforge.maths.integration.GaussRuleIntegrator; import hep.dataforge.maths.integration.UnivariateIntegrator; import hep.dataforge.plots.PlotFrame; import hep.dataforge.plots.data.PlottableXYFunction; +import hep.dataforge.utils.Misc; import hep.dataforge.values.NamedValueSet; import org.apache.commons.math3.analysis.BivariateFunction; import org.apache.commons.math3.analysis.UnivariateFunction; @@ -48,8 +49,10 @@ public class LossCalculator { private static final LossCalculator instance = new LossCalculator(); private static final UnivariateIntegrator integrator = new GaussRuleIntegrator(100); + private static final Map> lossProbCache = Misc.getLRUCache(100); private final Map cache = new HashMap<>(); + private LossCalculator() { cache.put(1, getSingleScatterFunction()); // cache.put(2, getDoubleScatterFunction()); @@ -281,36 +284,37 @@ public class LossCalculator { * рекурсивно вычисляем все вероятности, котрорые выше порога *

* дисер, стр.48 - * + *

* @param X * @return */ public List getLossProbabilities(double X) { + return lossProbCache.computeIfAbsent(X, (x) -> { + List res = new ArrayList<>(); + double prob; + if (x > 0) { + prob = 1 / x * (1 - Math.exp(-x)); + } else { + // если x ==0, то выживает только нулевой член, первый равен нулю + res.add(1d); + return res; + } + res.add(prob); - List res = new ArrayList<>(); - double prob; - if (X > 0) { - prob = 1 / X * (1 - Math.exp(-X)); - } else { - // если x ==0, то выживает только нулевой член, первый равен нулю - res.add(1d); - return res; - } - res.add(prob); - - while (prob > SCATTERING_PROBABILITY_THRESHOLD) { + while (prob > SCATTERING_PROBABILITY_THRESHOLD) { /* * prob(n) = prob(n-1)-1/n! * X^n * exp(-X); */ - double delta = Math.exp(-X); - for (int i = 1; i < res.size() + 1; i++) { - delta *= X / i; + double delta = Math.exp(-x); + for (int i = 1; i < res.size() + 1; i++) { + delta *= x / i; + } + prob -= delta / x; + res.add(prob); } - prob -= delta / X; - res.add(prob); - } - return res; + return res; + }); } public double getLossProbability(int order, double X) {