[no commit message]

This commit is contained in:
darksnake 2016-07-18 16:47:40 +03:00
parent bc05ae5006
commit 977331e8f8
6 changed files with 43 additions and 56 deletions

View File

@ -18,13 +18,14 @@ import inr.numass.actions.PileupSimulationAction
import hep.dataforge.grind.GrindMetaBuilder
import hep.dataforge.io.ColumnedDataWriter
File dataDir = new File("D:\\Work\\Numass\\data\\2016_04\\T2_data\\Fill_2_2\\set_6_e26d123e54010000")
if(!dataDir.exists()){
println "dataDir directory does not exist"
}
//File dataDir = new File("D:\\Work\\Numass\\data\\2016_04\\T2_data\\Fill_2_2\\set_7_b2a3433e54010000")
//File dataDir = new File("D:\\Work\\Numass\\data\\2016_04\\T2_data\\Fill_2_2\\set_6_e26d123e54010000")
//if(!dataDir.exists()){
// println "dataDir directory does not exist"
//}
//NumassData data = NumassDataLoader.fromLocalDir(null, dataDir)
NumassData data = NMFile.readFile(new File("D:\\Work\\Numass\\sterilie2013-2014\\dat\\2013\\SCAN06.DAT" ))
//println config
NumassData data = NumassDataLoader.fromLocalDir(null, dataDir)
Table t = new UnderflowCorrection().fitAllPoints(data, 500, 800, 20);
Table t = new UnderflowCorrection().fitAllPoints(data, 700, 1000,1800, 20);
ColumnedDataWriter.writeDataSet(System.out, t, "underflow parameters")

View File

@ -46,6 +46,7 @@ import inr.numass.models.EmpiricalLossSpectrum;
import inr.numass.models.ExperimentalVariableLossSpectrum;
import inr.numass.models.GaussSourceSpectrum;
import inr.numass.models.GunSpectrum;
import inr.numass.models.LossCalculator;
import inr.numass.models.ModularSpectrum;
import inr.numass.models.NBkgSpectrum;
import inr.numass.models.RangedNamedSetSpectrum;
@ -220,11 +221,18 @@ public class NumassPlugin extends BasicPlugin {
sp.setCaching(true);
}
//Adding trapping energy dependence
//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));
sp.setTrappingFunction((Ei, Ef) -> {
return 6.2e-5 * FastMath.exp(-(Ei - Ef) / 350d) + 1.97e-4 - 6.818e-9 * Ei;
});
switch (an.getString("trappingFunction", "default")) {
case "run2016":
sp.setTrappingFunction((Ei, Ef) -> {
return 6.2e-5 * FastMath.exp(-(Ei - Ef) / 350d) + 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);

View File

@ -50,12 +50,12 @@ import java.util.function.Function;
@ValueDef(name = "lowerWindowSlope", type = "NUMBER", def = "0", info = "Slope for the window lowerWindow bound")
@ValueDef(name = "upperWindow", type = "NUMBER", info = "Upper bound for window")
@ValueDef(name = "deadTime", type = "[NUMBER, STRING]", def = "0", info = "Dead time in s. Could be an expression.")
@ValueDef(name = "underflow", type = "BOOLEAN", def = "true",
info = "Enables calculation of detector threshold underflow using exponential shape of energy spectrum tail. "
+ "Not recomended to use with floating window.")
@ValueDef(name = "underflow.upperBorder", type = "NUMBER", def = "800", info = "Upper chanel for underflow calculation.")
@ValueDef(name = "underflow.threshold", type = "NUMBER", def = "17000", info = "The maximum U for undeflow calculation")
@ValueDef(name = "underflow.function", info = "An expression for underflow correction above threshold")
//@ValueDef(name = "underflow", type = "BOOLEAN", def = "true",
// info = "Enables calculation of detector threshold underflow using exponential shape of energy spectrum tail. "
// + "Not recomended to use with floating window.")
//@ValueDef(name = "underflow.upperBorder", type = "NUMBER", def = "800", info = "Upper chanel for underflow calculation.")
//@ValueDef(name = "underflow.threshold", type = "NUMBER", def = "17000", info = "The maximum U for undeflow calculation")
//@ValueDef(name = "underflow.function", info = "An expression for underflow correction above threshold")
@ValueDef(name = "correction",
info = "An expression to correct count tumber depending on potential ${U}, point length ${T} and point itself as ${point}")
public class PrepareDataAction extends OneToOneAction<NumassData, Table> {

View File

@ -17,7 +17,7 @@ import org.codehaus.groovy.control.customizers.ImportCustomizer;
*/
public class ExpressionUtils {
public static Double evaluate(String expression, Map<String, Object> binding) {
public static double evaluate(String expression, Map<String, Object> binding) {
Binding b = new Binding(binding);
// Add imports for script.
ImportCustomizer importCustomizer = new ImportCustomizer();
@ -28,6 +28,6 @@ public class ExpressionUtils {
configuration.addCompilationCustomizers(importCustomizer); // Create shell and execute script.
GroovyShell shell = new GroovyShell(b,configuration);
return (Double) shell.evaluate(expression);
return ((Number) shell.evaluate(expression)).doubleValue();
}
}

View File

@ -26,37 +26,6 @@ import java.io.FileNotFoundException;
import java.util.Locale;
import static java.util.Locale.setDefault;
import java.util.Scanner;
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;
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;
/**
*

View File

@ -28,7 +28,6 @@ public class UnderflowCorrection {
public double get(Reportable log, Meta meta, NMPoint point) {
if (point.getUset() >= meta.getDouble("underflow.threshold", 17000)) {
// log.report("Using underflow factor from formula: {}", meta.getString("underflow.function"));
if (meta.hasValue("underflow.function")) {
return TritiumUtils.evaluateExpression(point, meta.getString("underflow.function"));
} else {
@ -36,16 +35,13 @@ public class UnderflowCorrection {
}
} else {
try {
// log.report("Calculating underflow correction coefficient for point {}", point.getUset());
int xLow = meta.getInt("underflow.lowerBorder", meta.getInt("lowerWindow"));
int xHigh = meta.getInt("underflow.upperBorder", 800);
int binning = meta.getInt("underflow.binning", 20);
int upper = meta.getInt("upperWindow", RawNMPoint.MAX_CHANEL - 1);
long norm = point.getCountInWindow(xLow, upper);
double[] fitRes = getUnderflowExpParameters(point, xLow, xHigh, binning);
// log.report("Underflow interpolation function: {}*exp(c/{})", fitRes[0], fitRes[1]);
double correction = fitRes[0] * fitRes[1] * (Math.exp(xLow / fitRes[1]) - 1d) / norm + 1d;
// log.report("Underflow correction factor: {}", correction);
return correction;
} catch (Exception ex) {
log.reportError("Failed to calculate underflow parameters for point {} with message:", point.getUset(), ex.getMessage());
@ -62,6 +58,16 @@ public class UnderflowCorrection {
}
return builder.build();
}
public Table fitAllPoints(NumassData data, int xLow, int xHigh, int upper, int binning) {
ListTable.Builder builder = new ListTable.Builder("U", "amp", "expConst", "correction");
for (NMPoint point : data.getNMPoints()) {
long norm = point.getCountInWindow(xLow, upper);
double[] fitRes = getUnderflowExpParameters(point, xLow, xHigh, binning);
builder.row(point.getUset(), fitRes[0], fitRes[1], fitRes[0] * fitRes[1] * (Math.exp(xLow / fitRes[1]) - 1d) / norm + 1d);
}
return builder.build();
}
/**
* Calculate underflow exponent parameters using (xLow, xHigh) window for
@ -80,7 +86,10 @@ public class UnderflowCorrection {
List<WeightedObservedPoint> points = point.getMapWithBinning(binning, false)
.entrySet().stream()
.filter(entry -> entry.getKey() >= xLow && entry.getKey() <= xHigh)
.map(p -> new WeightedObservedPoint(1d / p.getValue(), p.getKey(), p.getValue() / binning))
.map(p -> new WeightedObservedPoint(
1d / p.getValue() * point.getLength() * point.getLength(), //weight
p.getKey(), // x
p.getValue() / binning / point.getLength())) //y
.collect(Collectors.toList());
SimpleCurveFitter fitter = SimpleCurveFitter.create(new ExponentFunction(), new double[]{1d, 200d});
return fitter.fit(points);