[no commit message]

This commit is contained in:
Alexander Nozik 2016-07-30 19:24:49 +03:00
parent bdd81ea455
commit 75e026d44f
19 changed files with 388 additions and 232 deletions

View File

@ -17,17 +17,17 @@ package inr.numass;
import hep.dataforge.actions.ActionUtils; import hep.dataforge.actions.ActionUtils;
import hep.dataforge.context.Context; import hep.dataforge.context.Context;
import hep.dataforge.context.GlobalContext;
import static hep.dataforge.context.GlobalContext.out; import static hep.dataforge.context.GlobalContext.out;
import hep.dataforge.data.FileDataFactory; import hep.dataforge.data.FileDataFactory;
import hep.dataforge.fitting.MINUITPlugin; import hep.dataforge.fitting.MINUITPlugin;
import hep.dataforge.io.IOManager; import hep.dataforge.io.IOManager;
import hep.dataforge.io.MetaFileReader; import hep.dataforge.io.MetaFileReader;
import hep.dataforge.meta.Meta; import hep.dataforge.meta.Meta;
import static inr.numass.NumassContext.printDescription; import static inr.numass.Numass.printDescription;
import java.io.File; import java.io.File;
import java.io.FileNotFoundException; import java.io.FileNotFoundException;
import java.util.Locale; import java.util.Locale;
import static java.util.Locale.setDefault;
import javax.swing.JFileChooser; import javax.swing.JFileChooser;
import javax.swing.JFrame; import javax.swing.JFrame;
import javax.swing.filechooser.FileFilter; import javax.swing.filechooser.FileFilter;
@ -41,20 +41,6 @@ import org.apache.commons.cli.ParseException;
import org.slf4j.Logger; import org.slf4j.Logger;
import org.slf4j.LoggerFactory; import org.slf4j.LoggerFactory;
import static java.util.Locale.setDefault; import static java.util.Locale.setDefault;
import static java.util.Locale.setDefault;
import static java.util.Locale.setDefault;
import static java.util.Locale.setDefault;
import static java.util.Locale.setDefault;
import static java.util.Locale.setDefault;
import static java.util.Locale.setDefault;
import static java.util.Locale.setDefault;
import static java.util.Locale.setDefault;
import static java.util.Locale.setDefault;
import static java.util.Locale.setDefault;
import static java.util.Locale.setDefault;
import static java.util.Locale.setDefault;
import static java.util.Locale.setDefault;
import static java.util.Locale.setDefault;
/** /**
* *
@ -64,13 +50,13 @@ public class Main {
public static void main(String[] args) throws Exception { public static void main(String[] args) throws Exception {
setDefault(Locale.US); setDefault(Locale.US);
NumassContext context = new NumassContext(); Context context = Numass.buildContext();
context.loadPlugin(new MINUITPlugin()); context.loadPlugin(new MINUITPlugin());
run(context, args); run(context, args);
} }
@SuppressWarnings("deprecation") @SuppressWarnings("deprecation")
public static void run(NumassContext context, String[] args) throws Exception { public static void run(Context context, String[] args) throws Exception {
Logger logger = LoggerFactory.getLogger("numass-main"); Logger logger = LoggerFactory.getLogger("numass-main");
Options options = prepareOptions(); Options options = prepareOptions();

View File

@ -30,27 +30,18 @@ import java.io.PrintWriter;
* *
* @author Darksnake * @author Darksnake
*/ */
public class NumassContext extends Context { public class Numass {
public NumassContext(Context parent, Meta config) { public static Context buildContext(Context parent, Meta meta) {
super(parent, "numass", config); Context numassContext = new Context(parent, "numass", meta);
init(); GlobalContext.registerContext(numassContext);
numassContext.loadPlugin("inr.numass:numass");
numassContext.setIO(new NumassIO());
return numassContext;
} }
public NumassContext(Context parent) { public static Context buildContext() {
super(parent, "numass"); return buildContext(GlobalContext.instance(), Meta.empty());
init();
}
public NumassContext() {
super("numass");
init();
}
private void init() {
GlobalContext.registerContext(this);
loadPlugin("inr.numass:numass");
setIO(new NumassIO());
} }
public static void printDescription(Context context, boolean allowANSI) throws DescriptorException { public static void printDescription(Context context, boolean allowANSI) throws DescriptorException {

View File

@ -48,7 +48,6 @@ import inr.numass.models.EmpiricalLossSpectrum;
import inr.numass.models.ExperimentalVariableLossSpectrum; import inr.numass.models.ExperimentalVariableLossSpectrum;
import inr.numass.models.GaussSourceSpectrum; import inr.numass.models.GaussSourceSpectrum;
import inr.numass.models.GunSpectrum; import inr.numass.models.GunSpectrum;
import inr.numass.models.LossCalculator;
import inr.numass.models.ModularSpectrum; import inr.numass.models.ModularSpectrum;
import inr.numass.models.NBkgSpectrum; import inr.numass.models.NBkgSpectrum;
import inr.numass.models.RangedNamedSetSpectrum; import inr.numass.models.RangedNamedSetSpectrum;
@ -98,11 +97,18 @@ public class NumassPlugin extends BasicPlugin {
} }
private void loadMath(MathPlugin math){ private void loadMath(MathPlugin math) {
math.registerBivariate("numass.trap2016", (Ei, Ef) -> {
return 3.92e-5 * FastMath.exp(-(Ei - Ef) / 300d) + 1.97e-4 - 6.818e-9 * Ei;
});
math.registerBivariate("numass.resolutionTail", meta -> {
double alpha = meta.getDouble("tailAlpha", 0);
double beta = meta.getDouble("tailBeta", 0);
return (double E, double U) -> 1 - (E - U) * (alpha + E / 1000d * beta) / 1000d;
});
} }
/** /**
* Load all numass model factories * Load all numass model factories
* *
@ -213,7 +219,7 @@ public class NumassPlugin extends BasicPlugin {
}); });
manager.addModel("sterile", (context, meta) -> { manager.addModel("sterile", (context, meta) -> {
SterileNeutrinoSpectrum sp = new SterileNeutrinoSpectrum(meta); SterileNeutrinoSpectrum sp = new SterileNeutrinoSpectrum(context, meta);
NBkgSpectrum spectrum = new NBkgSpectrum(sp); NBkgSpectrum spectrum = new NBkgSpectrum(sp);
return new XYModel(spectrum, getAdapter(meta)); return new XYModel(spectrum, getAdapter(meta));
@ -238,19 +244,11 @@ public class NumassPlugin extends BasicPlugin {
} }
//Adding trapping energy dependence //Adding trapping energy dependence
switch (meta.getString("trappingFunction", "default")) { if(meta.hasValue("transmission.trapping")){
case "run2016": BivariateFunction trap = MathPlugin.buildFrom(context).buildBivariateFunction(meta.getString("transmussion.trapping"));
sp.setTrappingFunction((Ei, Ef) -> { sp.setTrappingFunction(trap);
return 3.92e-5 * FastMath.exp(-(Ei - Ef) / 300d) + 1.97e-4 - 6.818e-9 * Ei;
});
break;
default:
//Intercept = 4.95745, B1 = -0.36879, B2 = 0.00827
sp.setTrappingFunction((Ei, Ef) -> LossCalculator.getTrapFunction().value(Ei, Ef) * (4.95745 - 0.36879 * Ei + 0.00827 * Ei * Ei));
} }
context.getReport().report("Using folowing trapping formula: {}",
"6.2e-5 * FastMath.exp(-(Ei - Ef) / 350d) + 1.97e-4 - 6.818e-9 * Ei");
NBkgSpectrum spectrum = new NBkgSpectrum(sp); NBkgSpectrum spectrum = new NBkgSpectrum(sp);
return new XYModel(spectrum, getAdapter(meta)); return new XYModel(spectrum, getAdapter(meta));

View File

@ -23,7 +23,7 @@ public class WorkspaceTest {
String storagepath = "D:\\Work\\Numass\\data\\"; String storagepath = "D:\\Work\\Numass\\data\\";
Workspace workspace = BasicWorkspace.builder() Workspace workspace = BasicWorkspace.builder()
.setContext(new NumassContext()) .setContext(Numass.buildContext())
.loadData("", new StorageDataFactory(), new MetaBuilder("storage").putValue("path", storagepath)) .loadData("", new StorageDataFactory(), new MetaBuilder("storage").putValue("path", storagepath))
.build(); .build();

View File

@ -41,7 +41,7 @@ import hep.dataforge.tables.Table;
import hep.dataforge.tables.XYAdapter; import hep.dataforge.tables.XYAdapter;
import hep.dataforge.values.NamedValueSet; import hep.dataforge.values.NamedValueSet;
import inr.numass.NumassIntegrator; import inr.numass.NumassIntegrator;
import inr.numass.NumassContext; import inr.numass.Numass;
import inr.numass.models.ExperimentalVariableLossSpectrum; import inr.numass.models.ExperimentalVariableLossSpectrum;
import inr.numass.models.LossCalculator; import inr.numass.models.LossCalculator;
import java.io.OutputStreamWriter; import java.io.OutputStreamWriter;

View File

@ -8,7 +8,7 @@ package inr.numass.models;
import hep.dataforge.fitting.parametric.ParametricFunction; import hep.dataforge.fitting.parametric.ParametricFunction;
import hep.dataforge.values.NamedValueSet; import hep.dataforge.values.NamedValueSet;
import inr.numass.NumassIntegrator; import inr.numass.NumassIntegrator;
import inr.numass.NumassContext; import inr.numass.Numass;
import inr.numass.utils.TritiumUtils; import inr.numass.utils.TritiumUtils;
import org.apache.commons.math3.analysis.UnivariateFunction; import org.apache.commons.math3.analysis.UnivariateFunction;

View File

@ -326,7 +326,7 @@ public class LossCalculator {
} }
public synchronized double getLossValue(int order, double Ei, double Ef) { public synchronized double getLossValue(int order, double Ei, double Ef) {
if (Ei < Ef) { if (Ei - Ef < 5d) {
return 0; return 0;
} else if (Ei - Ef >= getMargin(order)) { } else if (Ei - Ef >= getMargin(order)) {
return 0; return 0;
@ -377,7 +377,7 @@ public class LossCalculator {
return 0; return 0;
} }
}; };
return integrator.integrate(integrand, 0d, margin); return integrator.integrate(integrand, 5d, margin);
}; };
return FunctionCaching.cacheUnivariateFunction(res, 0, margin, 200); return FunctionCaching.cacheUnivariateFunction(res, 0, margin, 200);

View File

@ -16,7 +16,6 @@
package inr.numass.models; package inr.numass.models;
import inr.numass.NumassIntegrator; import inr.numass.NumassIntegrator;
import inr.numass.NumassContext;
import org.apache.commons.math3.analysis.BivariateFunction; import org.apache.commons.math3.analysis.BivariateFunction;
import org.apache.commons.math3.analysis.UnivariateFunction; import org.apache.commons.math3.analysis.UnivariateFunction;

View File

@ -19,7 +19,7 @@ import hep.dataforge.fitting.parametric.AbstractParametricFunction;
import hep.dataforge.fitting.parametric.ParametricFunction; import hep.dataforge.fitting.parametric.ParametricFunction;
import hep.dataforge.values.NamedValueSet; import hep.dataforge.values.NamedValueSet;
import inr.numass.NumassIntegrator; import inr.numass.NumassIntegrator;
import inr.numass.NumassContext; import inr.numass.Numass;
import org.apache.commons.math3.analysis.BivariateFunction; import org.apache.commons.math3.analysis.BivariateFunction;
import org.apache.commons.math3.analysis.UnivariateFunction; import org.apache.commons.math3.analysis.UnivariateFunction;

View File

@ -22,7 +22,7 @@ import hep.dataforge.maths.integration.UnivariateIntegrator;
import hep.dataforge.values.NamedValueSet; import hep.dataforge.values.NamedValueSet;
import hep.dataforge.values.ValueProvider; import hep.dataforge.values.ValueProvider;
import inr.numass.NumassIntegrator; import inr.numass.NumassIntegrator;
import inr.numass.NumassContext; import inr.numass.Numass;
import java.util.List; import java.util.List;
import org.apache.commons.math3.analysis.BivariateFunction; import org.apache.commons.math3.analysis.BivariateFunction;
import org.apache.commons.math3.analysis.UnivariateFunction; import org.apache.commons.math3.analysis.UnivariateFunction;

View File

@ -9,7 +9,6 @@ import hep.dataforge.fitting.parametric.AbstractParametricBiFunction;
import hep.dataforge.meta.Meta; import hep.dataforge.meta.Meta;
import hep.dataforge.values.NamedValueSet; import hep.dataforge.values.NamedValueSet;
import inr.numass.models.ResolutionFunction; import inr.numass.models.ResolutionFunction;
import static java.lang.Double.isNaN;
import static java.lang.Math.sqrt; import static java.lang.Math.sqrt;
import org.apache.commons.math3.analysis.BivariateFunction; import org.apache.commons.math3.analysis.BivariateFunction;

View File

@ -5,7 +5,9 @@
*/ */
package inr.numass.models.sterile; package inr.numass.models.sterile;
import hep.dataforge.context.Context;
import hep.dataforge.fitting.parametric.AbstractParametricBiFunction; import hep.dataforge.fitting.parametric.AbstractParametricBiFunction;
import hep.dataforge.maths.MathPlugin;
import hep.dataforge.meta.Meta; import hep.dataforge.meta.Meta;
import hep.dataforge.values.NamedValueSet; import hep.dataforge.values.NamedValueSet;
import inr.numass.models.LossCalculator; import inr.numass.models.LossCalculator;
@ -22,21 +24,18 @@ public class NumassTransmission extends AbstractParametricBiFunction {
private final LossCalculator calculator; private final LossCalculator calculator;
private final BivariateFunction trapFunc; private final BivariateFunction trapFunc;
// private BicubicInterpolatingFunction cache; public NumassTransmission(Context context, Meta meta) {
// private double cachedX;
// private Meta cacheMeta;
public NumassTransmission(Meta meta) {
super(list); super(list);
this.calculator = LossCalculator.instance(); this.calculator = LossCalculator.instance();
trapFunc = getTrapFunction(meta.getNodeOrEmpty("trapping")); if (meta.hasValue("trapping")) {
// if (meta.hasNode("cache")) { trapFunc = MathPlugin.buildFrom(context).buildBivariateFunction(meta.getString("trapping"));
// cacheMeta = meta.getNode("cache"); } else {
// } trapFunc = getTrapFunction(context, meta.getNodeOrEmpty("trapping"));
}
} }
private BivariateFunction getTrapFunction(Meta meta) { private BivariateFunction getTrapFunction(Context context, Meta meta) {
return LossCalculator.getTrapFunction(); return MathPlugin.buildFrom(context).buildBivariateFunction(meta);
} }
@Override @Override
@ -65,18 +64,6 @@ public class NumassTransmission extends AbstractParametricBiFunction {
return LossCalculator.instance().getLossProbability(0, getX(eIn, set)); return LossCalculator.instance().getLossProbability(0, getX(eIn, set));
} }
// private synchronized void setupCache(double X) {
// if (this.cachedX != X) {
// double cacheLo = cacheMeta.getDouble("lo", 14000);
// double cacheHi = cacheMeta.getDouble("hi", 18575);
// int numPoints = cacheMeta.getInt("numPoints", 1000);
// double[] eIns = GridCalculator.getUniformUnivariateGrid(cacheLo, cacheHi, numPoints);
// double[] eOuts = GridCalculator.getUniformUnivariateGrid(cacheLo, cacheHi, numPoints);
// double[][] vals = MathUtils.calculateFunction(calculator.getTotalLossBivariateFunction(X), eOuts, eOuts);
// this.cachedX = X;
// this.cache = new BicubicInterpolator().interpolate(eIns, eOuts, vals);
// }
// }
@Override @Override
public double value(double eIn, double eOut, NamedValueSet set) { public double value(double eIn, double eOut, NamedValueSet set) {
//calculate X taking into account its energy dependence //calculate X taking into account its energy dependence

View File

@ -10,6 +10,7 @@ import hep.dataforge.context.GlobalContext;
import hep.dataforge.description.NodeDef; import hep.dataforge.description.NodeDef;
import hep.dataforge.description.ValueDef; import hep.dataforge.description.ValueDef;
import hep.dataforge.exceptions.NotDefinedException; import hep.dataforge.exceptions.NotDefinedException;
import hep.dataforge.fitting.parametric.AbstractParametricBiFunction;
import hep.dataforge.fitting.parametric.AbstractParametricFunction; import hep.dataforge.fitting.parametric.AbstractParametricFunction;
import hep.dataforge.fitting.parametric.ParametricBiFunction; import hep.dataforge.fitting.parametric.ParametricBiFunction;
import hep.dataforge.maths.integration.UnivariateIntegrator; import hep.dataforge.maths.integration.UnivariateIntegrator;
@ -17,14 +18,7 @@ import hep.dataforge.meta.Meta;
import hep.dataforge.values.NamedValueSet; import hep.dataforge.values.NamedValueSet;
import inr.numass.NumassIntegrator; import inr.numass.NumassIntegrator;
import inr.numass.models.FSS; import inr.numass.models.FSS;
import java.util.stream.DoubleStream;
import org.apache.commons.math3.analysis.BivariateFunction;
import org.apache.commons.math3.analysis.UnivariateFunction; import org.apache.commons.math3.analysis.UnivariateFunction;
import org.apache.commons.math3.distribution.EnumeratedRealDistribution;
import org.apache.commons.math3.distribution.RealDistribution;
import org.apache.commons.math3.random.JDKRandomGenerator;
import org.apache.commons.math3.random.RandomGenerator;
import org.apache.commons.math3.random.SynchronizedRandomGenerator;
/** /**
* Compact all-in-one model for sterile neutrino spectrum * Compact all-in-one model for sterile neutrino spectrum
@ -39,8 +33,8 @@ public class SterileNeutrinoSpectrum extends AbstractParametricFunction {
private static final String[] list = {"X", "trap", "E0", "mnu2", "msterile2", "U2"}; private static final String[] list = {"X", "trap", "E0", "mnu2", "msterile2", "U2"};
private final RandomGenerator rnd; // private final RandomGenerator rnd;
private RealDistribution fssDistribution; // private RealDistribution fssDistribution;
private FSS fss; private FSS fss;
/** /**
@ -57,21 +51,27 @@ public class SterileNeutrinoSpectrum extends AbstractParametricFunction {
*/ */
private final ParametricBiFunction resolution; private final ParametricBiFunction resolution;
private boolean useMC; /**
* auxiliary function for trans-res convolution
*/
private final ParametricBiFunction transRes;
// private boolean useMC;
private boolean fast; private boolean fast;
public SterileNeutrinoSpectrum(Context context, Meta configuration) { public SterileNeutrinoSpectrum(Context context, Meta configuration) {
super(list); super(list);
rnd = new SynchronizedRandomGenerator(new JDKRandomGenerator()); // rnd = new SynchronizedRandomGenerator(new JDKRandomGenerator());
if (configuration.hasValue("fssFile")) { if (configuration.hasValue("fssFile")) {
fss = new FSS(context.io().getFile(configuration.getString("fssFile"))); fss = new FSS(context.io().getFile(configuration.getString("fssFile")));
fssDistribution = new EnumeratedRealDistribution(rnd, fss.getEs(), fss.getPs()); // fssDistribution = new EnumeratedRealDistribution(rnd, fss.getEs(), fss.getPs());
} }
transmission = new NumassTransmission(configuration.getNodeOrEmpty("transmission")); transmission = new NumassTransmission(context, configuration.getNodeOrEmpty("transmission"));
resolution = new NumassResolution(configuration.getNode("resolution", Meta.empty())); resolution = new NumassResolution(configuration.getNode("resolution", Meta.empty()));
this.useMC = configuration.getBoolean("useMC", false); // this.useMC = configuration.getBoolean("useMC", false);
this.fast = configuration.getBoolean("fast", true); this.fast = configuration.getBoolean("fast", true);
transRes = new TransRes();
} }
public SterileNeutrinoSpectrum(Meta configuration) { public SterileNeutrinoSpectrum(Meta configuration) {
@ -89,10 +89,10 @@ public class SterileNeutrinoSpectrum extends AbstractParametricFunction {
case "msterile2": case "msterile2":
case "mnu2": case "mnu2":
case "E0": case "E0":
return integrate(u, source.derivative(parName), transmission, resolution, set); return integrate(u, source.derivative(parName), transRes, set);
case "X": case "X":
case "trap": case "trap":
return integrate(u, source, transmission, resolution.derivative(parName), set); return integrate(u, source, transRes.derivative(parName), set);
default: default:
throw new NotDefinedException(); throw new NotDefinedException();
} }
@ -100,39 +100,94 @@ public class SterileNeutrinoSpectrum extends AbstractParametricFunction {
@Override @Override
public double value(double u, NamedValueSet set) { public double value(double u, NamedValueSet set) {
if (useDirect()) { return integrate(u, source, transRes, set);
return integrateDirect(u, source, transmission, resolution, set);
} else {
return integrate(u, source, transmission, resolution, set);
}
}
private int numCalls(double u) {
return 100000;
}
private boolean useDirect() {
return !useMC;
} }
// private int numCalls(double u) {
// return 100000;
// }
//
// private boolean useDirect() {
// return !useMC;
// }
@Override @Override
public boolean providesDeriv(String name) { public boolean providesDeriv(String name) {
return source.providesDeriv(name) && transmission.providesDeriv(name) && resolution.providesDeriv(name); return source.providesDeriv(name) && transmission.providesDeriv(name) && resolution.providesDeriv(name);
} }
// /**
// * Random E generator
// *
// * @param a
// * @param b
// * @return
// */
// private double rndE(double a, double b) {
// return rnd.nextDouble() * (b - a) + a;
// }
//
// private double integrate(
// double u,
// ParametricBiFunction sourceFunction,
// ParametricBiFunction transmissionFunction,
// ParametricBiFunction resolutionFunction,
// NamedValueSet set) {
// if (useDirect()) {
// return integrateDirect(u, sourceFunction, transmissionFunction, resolutionFunction, set);
// } else {
// return integrateRandom(u, sourceFunction, transmissionFunction, resolutionFunction, set);
// }
// }
// /**
// * Monte-Carlo integration of spectrum
// *
// * @param u
// * @param sourceFunction
// * @param transmissionFunction
// * @param resolutionFunction
// * @param set
// * @return
// */
// private double integrateRandom(
// double u,
// ParametricBiFunction sourceFunction,
// ParametricBiFunction transmissionFunction,
// ParametricBiFunction resolutionFunction,
// NamedValueSet set) {
//
// int num = numCalls(u);
// double eMax = set.getDouble("E0") + 5d;
// if (u > eMax) {
// return 0;
// }
//
// double sum = DoubleStream.generate(() -> {
// // generate final state
// double fs;
// if (fssDistribution != null) {
// fs = fssDistribution.sample();
// } else {
// fs = 0;
// }
//
// double eIn = rndE(u, eMax);
//
// double eOut = rndE(u, eIn);
//
// double res = sourceFunction.value(fs, eIn, set)
// * transmissionFunction.value(eIn, eOut, set)
// * resolutionFunction.value(eOut, u, set);
//
// if (Double.isNaN(res)) {
// throw new Error();
// }
// return res;
// }).parallel().limit(num).sum();
// //triangle surface
// return Math.pow(eMax - u, 2d) / 2d * sum / num;
// }
/** /**
* Random E generator * Direct Gauss-Legandre integration
*
* @param a
* @param b
* @return
*/
private double rndE(double a, double b) {
return rnd.nextDouble() * (b - a) + a;
}
/**
* Monte-Carlo integration of spectrum
* *
* @param u * @param u
* @param sourceFunction * @param sourceFunction
@ -144,57 +199,7 @@ public class SterileNeutrinoSpectrum extends AbstractParametricFunction {
private double integrate( private double integrate(
double u, double u,
ParametricBiFunction sourceFunction, ParametricBiFunction sourceFunction,
ParametricBiFunction transmissionFunction, ParametricBiFunction transResFunction,
ParametricBiFunction resolutionFunction,
NamedValueSet set) {
int num = numCalls(u);
double eMax = set.getDouble("E0") + 5d;
if (u > eMax) {
return 0;
}
double sum = DoubleStream.generate(() -> {
// generate final state
double fs;
if (fssDistribution != null) {
fs = fssDistribution.sample();
} else {
fs = 0;
}
double eIn = rndE(u, eMax);
double eOut = rndE(u, eIn);
double res = sourceFunction.value(fs, eIn, set)
* transmissionFunction.value(eIn, eOut, set)
* resolutionFunction.value(eOut, u, set);
if (Double.isNaN(res)) {
throw new Error();
}
return res;
}).parallel().limit(num).sum();
//triangle surface
return Math.pow(eMax - u, 2d) / 2d * sum / num;
}
/**
* Direct Gauss-Legandre integration
*
* @param u
* @param sourceFunction
* @param transmissionFunction
* @param resolutionFunction
* @param set
* @return
*/
private double integrateDirect(
double u,
ParametricBiFunction sourceFunction,
ParametricBiFunction transmissionFunction,
ParametricBiFunction resolutionFunction,
NamedValueSet set) { NamedValueSet set) {
double eMax = set.getDouble("E0") + 5d; double eMax = set.getDouble("E0") + 5d;
@ -216,14 +221,55 @@ public class SterileNeutrinoSpectrum extends AbstractParametricFunction {
fsSource = (eIn) -> sourceFunction.value(0, eIn, set); fsSource = (eIn) -> sourceFunction.value(0, eIn, set);
} }
BivariateFunction transRes = (eIn, usp) -> { UnivariateIntegrator integrator;
UnivariateFunction integrand = (eOut) -> transmissionFunction.value(eIn, eOut, set) * resolutionFunction.value(eOut, usp, set); if (fast && eMax - u < 500) {
integrator = NumassIntegrator.getFastInterator();
} else {
integrator = NumassIntegrator.getDefaultIntegrator();
}
return integrator.integrate(eIn -> fsSource.value(eIn) * transResFunction.value(eIn, u, set), u, eMax);
}
private class TransRes extends AbstractParametricBiFunction {
public TransRes() {
super(new String[]{"X", "trap"});
}
@Override
public boolean providesDeriv(String name) {
return true;
}
@Override
public double derivValue(String parName, double eIn, double u, NamedValueSet set) {
switch (parName) {
case "X":
//TODO implement p0 derivative
throw new NotDefinedException();
case "trap":
return lossRes(transmission.derivative(parName), eIn, u, set);
default:
return super.derivValue(parName, eIn, u, set);
}
}
@Override
public double value(double eIn, double u, NamedValueSet set) {
double p0 = NumassTransmission.p0(eIn, set);
return p0 * resolution.value(eIn, u, set) + lossRes(transmission, eIn, u, set);
}
private double lossRes(ParametricBiFunction transFunc, double eIn, double u, NamedValueSet set) {
UnivariateFunction integrand = (eOut) -> transFunc.value(eIn, eOut, set) * resolution.value(eOut, u, set);
double border = u + 30; double border = u + 30;
double firstPart = NumassIntegrator.getFastInterator().integrate(integrand, usp, Math.min(eIn, border)); double firstPart = NumassIntegrator.getFastInterator().integrate(integrand, u, Math.min(eIn, border));
double secondPart; double secondPart;
if (eIn > border) { if (eIn > border) {
if(fast){ if (fast) {
secondPart = NumassIntegrator.getFastInterator().integrate(integrand, border, eIn); secondPart = NumassIntegrator.getFastInterator().integrate(integrand, border, eIn);
} else { } else {
secondPart = NumassIntegrator.getDefaultIntegrator().integrate(integrand, border, eIn); secondPart = NumassIntegrator.getDefaultIntegrator().integrate(integrand, border, eIn);
@ -232,21 +278,8 @@ public class SterileNeutrinoSpectrum extends AbstractParametricFunction {
secondPart = 0; secondPart = 0;
} }
return firstPart + secondPart; return firstPart + secondPart;
};
UnivariateFunction integrand = (eIn) -> {
double p0 = NumassTransmission.p0(eIn, set);
return fsSource.value(eIn) * (p0 * resolutionFunction.value(eIn, u, set) + transRes.value(eIn, u));
};
UnivariateIntegrator integrator;
if (fast && eMax - u < 500) {
integrator = NumassIntegrator.getFastInterator();
} else {
integrator = NumassIntegrator.getDefaultIntegrator();
} }
return integrator.integrate(integrand, u, eMax);
} }
} }

View File

@ -0,0 +1,102 @@
/*
* To change this license header, choose License Headers in Project Properties.
* To change this template file, choose Tools | Templates
* and open the template in the editor.
*/
package inr.numass.models.sterile;
import hep.dataforge.context.Context;
import hep.dataforge.fitting.ParamSet;
import hep.dataforge.fitting.parametric.ParametricFunction;
import hep.dataforge.maths.MathPlugin;
import hep.dataforge.meta.Meta;
import hep.dataforge.meta.MetaBuilder;
import inr.numass.Numass;
import inr.numass.models.BetaSpectrum;
import inr.numass.models.ModularSpectrum;
import inr.numass.models.NBkgSpectrum;
import inr.numass.models.RangedNamedSetSpectrum;
import inr.numass.models.ResolutionFunction;
import org.apache.commons.math3.analysis.BivariateFunction;
/**
*
* @author Alexander Nozik
*/
public class TestModels {
/**
* @param args the command line arguments
*/
public static void main(String[] args) {
Context context = Numass.buildContext();
/*
<model modelName="sterile" fssFile="FS.txt">
<resolution width = "1.22e-4" tailAlpha="1e-2"/>
<transmission trapping="numass.trap2016"/>
</model>
*/
Meta meta = new MetaBuilder("model")
.putNode(new MetaBuilder("resolution")
.putValue("width", 1.22e-4)
.putValue("tailAlpha", 1e-2)
)
.putNode(new MetaBuilder("transmission")
.putValue("trapping", "numass.trap2016")
);
ParametricFunction oldFunc = oldModel(context, meta);
ParametricFunction newFunc = newModel(context, meta);
ParamSet allPars = new ParamSet()
.setPar("N", 7e+05, 1.8e+03, 0d, Double.POSITIVE_INFINITY)
.setPar("bkg", 1d, 0.050)
.setPar("E0", 18575d, 1.4)
.setPar("mnu2", 0d, 1d)
.setPar("msterile2", 1000d * 1000d, 0)
.setPar("U2", 0.0, 1e-4, -1d, 1d)
.setPar("X", 0.04, 0.01, 0d, Double.POSITIVE_INFINITY)
.setPar("trap", 1, 0.01, 0d, Double.POSITIVE_INFINITY);
for (double u = 14000; u < 18600; u += 100) {
// double oldVal = oldFunc.value(u, allPars);
// double newVal = newFunc.value(u, allPars);
double oldVal = oldFunc.derivValue("trap", u, allPars);
double newVal = newFunc.derivValue("trap", u, allPars);
System.out.printf("%f\t%g\t%g\t%g%n", u, oldVal, newVal, 1d - oldVal / newVal);
}
}
private static ParametricFunction oldModel(Context context, Meta meta) {
double A = meta.getDouble("resolution", meta.getDouble("resolution.width", 8.3e-5));//8.3e-5
double from = meta.getDouble("from", 13900d);
double to = meta.getDouble("to", 18700d);
context.getReport().report("Setting up tritium model with real transmission function");
BivariateFunction resolutionTail;
if (meta.hasValue("resolution.tailAlpha")) {
resolutionTail = ResolutionFunction.getAngledTail(meta.getDouble("resolution.tailAlpha"), meta.getDouble("resolution.tailBeta", 0));
} else {
resolutionTail = ResolutionFunction.getRealTail();
}
//RangedNamedSetSpectrum beta = new BetaSpectrum(context.io().getFile("FS.txt"));
RangedNamedSetSpectrum beta = new BetaSpectrum();
ModularSpectrum sp = new ModularSpectrum(beta, new ResolutionFunction(A, resolutionTail), from, to);
if (meta.getBoolean("caching", false)) {
context.getReport().report("Caching turned on");
sp.setCaching(true);
}
//Adding trapping energy dependence
if (meta.hasValue("transmission.trapping")) {
BivariateFunction trap = MathPlugin.buildFrom(context).buildBivariateFunction(meta.getString("transmission.trapping"));
sp.setTrappingFunction(trap);
}
return new NBkgSpectrum(sp);
}
private static ParametricFunction newModel(Context context, Meta meta) {
SterileNeutrinoSpectrum sp = new SterileNeutrinoSpectrum(context, meta);
return new NBkgSpectrum(sp);
}
}

View File

@ -30,7 +30,6 @@ import hep.dataforge.meta.MetaBuilder;
import hep.dataforge.plots.PlotFrame; import hep.dataforge.plots.PlotFrame;
import hep.dataforge.plots.PlotHolder; import hep.dataforge.plots.PlotHolder;
import hep.dataforge.plots.PlotsPlugin; import hep.dataforge.plots.PlotsPlugin;
import hep.dataforge.utils.MetaFactory;
import hep.dataforge.values.Value; import hep.dataforge.values.Value;
import inr.numass.NumassIO; import inr.numass.NumassIO;
import inr.numass.NumassProperties; import inr.numass.NumassProperties;
@ -60,6 +59,7 @@ import javafx.scene.control.TitledPane;
import javafx.scene.control.ToggleButton; import javafx.scene.control.ToggleButton;
import javafx.stage.FileChooser; import javafx.stage.FileChooser;
import org.controlsfx.control.StatusBar; import org.controlsfx.control.StatusBar;
import hep.dataforge.utils.ContextMetaFactory;
/** /**
* FXML Controller class * FXML Controller class
@ -69,7 +69,7 @@ import org.controlsfx.control.StatusBar;
public class NumassWorkbenchController implements Initializable, StagePaneHolder, ActionStateListener, PlotHolder { public class NumassWorkbenchController implements Initializable, StagePaneHolder, ActionStateListener, PlotHolder {
Context parentContext; Context parentContext;
MetaFactory<Context> contextFactory; ContextMetaFactory<Context> contextFactory;
List<MetaEditor> actionEditors = new ArrayList<>(); List<MetaEditor> actionEditors = new ArrayList<>();
MetaEditor dataEditor; MetaEditor dataEditor;
@ -402,7 +402,7 @@ public class NumassWorkbenchController implements Initializable, StagePaneHolder
this.parentContext = parentContext; this.parentContext = parentContext;
} }
public void setContextFactory(MetaFactory<Context> contextFactory) { public void setContextFactory(ContextMetaFactory<Context> contextFactory) {
this.contextFactory = contextFactory; this.contextFactory = contextFactory;
} }

View File

@ -6,7 +6,7 @@
package inr.numass.workbench; package inr.numass.workbench;
import hep.dataforge.context.GlobalContext; import hep.dataforge.context.GlobalContext;
import inr.numass.NumassContext; import inr.numass.Numass;
import java.io.IOException; import java.io.IOException;
import java.text.ParseException; import java.text.ParseException;
import javafx.application.Application; import javafx.application.Application;
@ -30,7 +30,7 @@ public class Workbench extends Application {
Scene scene = new Scene(parent, 800, 600); Scene scene = new Scene(parent, 800, 600);
NumassWorkbenchController controller = loader.getController(); NumassWorkbenchController controller = loader.getController();
controller.setContextFactory(NumassContext::new); controller.setContextFactory(Numass::buildContext);
primaryStage.setTitle("Numass workbench"); primaryStage.setTitle("Numass workbench");
primaryStage.setScene(scene); primaryStage.setScene(scene);

View File

@ -0,0 +1,60 @@
/*
* To change this license header, choose License Headers in Project Properties.
* To change this template file, choose Tools | Templates
* and open the template in the editor.
*/
package inr.numass;
import hep.dataforge.context.Context;
import org.apache.commons.cli.CommandLine;
import org.junit.Test;
import static org.junit.Assert.*;
/**
*
* @author Alexander Nozik
*/
public class MainTest {
public MainTest() {
}
/**
* Test of main method, of class Main.
*/
@Test
public void testMain() throws Exception {
System.out.println("main");
String[] args = null;
Main.main(args);
// TODO review the generated test code and remove the default call to fail.
fail("The test case is a prototype.");
}
/**
* Test of run method, of class Main.
*/
@Test
public void testRun() throws Exception {
System.out.println("run");
Context context = null;
String[] args = null;
Main.run(context, args);
// TODO review the generated test code and remove the default call to fail.
fail("The test case is a prototype.");
}
/**
* Test of applyCLItoContext method, of class Main.
*/
@Test
public void testApplyCLItoContext() throws Exception {
System.out.println("applyCLItoContext");
CommandLine line = null;
Context context = null;
Main.applyCLItoContext(line, context);
// TODO review the generated test code and remove the default call to fail.
fail("The test case is a prototype.");
}
}

View File

@ -1,32 +0,0 @@
/*
* To change this license header, choose License Headers in Project Properties.
* To change this template file, choose Tools | Templates
* and open the template in the editor.
*/
package inr.numass;
import inr.numass.storage.SetDirectionUtility;
import org.junit.Test;
/**
*
* @author Alexander Nozik
*/
public class NumassContextTest {
public NumassContextTest() {
}
/**
* Test of close method, of class NumassContext.
*/
@Test
public void testClose() throws Exception {
System.out.println("close");
NumassContext instance = new NumassContext();
instance.close();
assert SetDirectionUtility.cacheFile(instance).exists();
}
}

View File

@ -0,0 +1,33 @@
/*
* To change this license header, choose License Headers in Project Properties.
* To change this template file, choose Tools | Templates
* and open the template in the editor.
*/
package inr.numass;
import hep.dataforge.context.Context;
import hep.dataforge.maths.MathPlugin;
import hep.dataforge.plots.PlotsPlugin;
import org.junit.Test;
/**
*
* @author Alexander Nozik
*/
public class NumassTest {
public NumassTest() {
}
/**
* Test of buildContext method, of class Numass.
*/
@Test
public void testBuildContext() {
Context context = Numass.buildContext();
MathPlugin.buildFrom(context);
PlotsPlugin.buildFrom(context);
}
}