diff --git a/numass-core/src/main/java/inr/numass/data/analyzers/AbstractAnalyzer.java b/numass-core/src/main/java/inr/numass/data/analyzers/AbstractAnalyzer.java index b3cebb37..d83d7af0 100644 --- a/numass-core/src/main/java/inr/numass/data/analyzers/AbstractAnalyzer.java +++ b/numass-core/src/main/java/inr/numass/data/analyzers/AbstractAnalyzer.java @@ -19,9 +19,10 @@ import static inr.numass.data.api.NumassPoint.HV_KEY; */ public abstract class AbstractAnalyzer implements NumassAnalyzer { public static String WINDOW_KEY = "window"; + public static String TIME_KEY = "timestamp"; - public static String[] NAME_LIST = {LENGTH_KEY, COUNT_KEY, COUNT_RATE_KEY, COUNT_RATE_ERROR_KEY, WINDOW_KEY, "timestamp"}; - public static String[] NAME_LIST_WITH_HV = {HV_KEY, LENGTH_KEY, COUNT_KEY, COUNT_RATE_KEY, COUNT_RATE_ERROR_KEY, WINDOW_KEY, "timestamp"}; + public static String[] NAME_LIST = {LENGTH_KEY, COUNT_KEY, COUNT_RATE_KEY, COUNT_RATE_ERROR_KEY, WINDOW_KEY, TIME_KEY}; + public static String[] NAME_LIST_WITH_HV = {HV_KEY, LENGTH_KEY, COUNT_KEY, COUNT_RATE_KEY, COUNT_RATE_ERROR_KEY, WINDOW_KEY, TIME_KEY}; @Nullable private final SignalProcessor processor; @@ -64,10 +65,13 @@ public abstract class AbstractAnalyzer implements NumassAnalyzer { } } - - @Override - public Table analyze(NumassSet set, Meta config) { - TableFormat format = new TableFormatBuilder() + /** + * Get table format for summary table + * @param config + * @return + */ + protected TableFormat getTableFormat(Meta config){ + return new TableFormatBuilder() .addNumber(HV_KEY, X_VALUE_KEY) .addNumber(LENGTH_KEY) .addNumber(COUNT_KEY) @@ -76,6 +80,12 @@ public abstract class AbstractAnalyzer implements NumassAnalyzer { .addColumn(WINDOW_KEY) .addTime() .build(); + } + + + @Override + public Table analyze(NumassSet set, Meta config) { + TableFormat format = getTableFormat(config); return new ListTable.Builder(format) .rows(set.getPoints().map(point -> analyze(point, config))) diff --git a/numass-core/src/main/java/inr/numass/data/analyzers/SimpleAnalyzer.java b/numass-core/src/main/java/inr/numass/data/analyzers/SimpleAnalyzer.java index 10780815..45cacfda 100644 --- a/numass-core/src/main/java/inr/numass/data/analyzers/SimpleAnalyzer.java +++ b/numass-core/src/main/java/inr/numass/data/analyzers/SimpleAnalyzer.java @@ -27,13 +27,14 @@ public class SimpleAnalyzer extends AbstractAnalyzer { int loChannel = config.getInt("window.lo", 0); int upChannel = config.getInt("window.up", Integer.MAX_VALUE); long count = getEvents(block, config).count(); - double countRate = (double) count / block.getLength().toMillis() * 1000; - double countRateError = Math.sqrt((double) count) / block.getLength().toMillis() * 1000; + double length = (double) block.getLength().toNanos()/1e9; + double countRate = (double) count / length; + double countRateError = Math.sqrt((double) count) / length; if (block instanceof NumassPoint) { return ValueMap.of(NAME_LIST_WITH_HV, ((NumassPoint) block).getVoltage(), - block.getLength().toNanos(), + length, count, countRate, countRateError, @@ -41,7 +42,7 @@ public class SimpleAnalyzer extends AbstractAnalyzer { block.getStartTime()); } else { return ValueMap.of(NAME_LIST, - block.getLength().toNanos(), + length, count, countRate, countRateError, diff --git a/numass-core/src/main/java/inr/numass/data/analyzers/SmartAnalyzer.java b/numass-core/src/main/java/inr/numass/data/analyzers/SmartAnalyzer.java index 6c61e4a9..ebf64bd3 100644 --- a/numass-core/src/main/java/inr/numass/data/analyzers/SmartAnalyzer.java +++ b/numass-core/src/main/java/inr/numass/data/analyzers/SmartAnalyzer.java @@ -1,11 +1,18 @@ package inr.numass.data.analyzers; +import hep.dataforge.description.ValueDef; import hep.dataforge.meta.Meta; +import hep.dataforge.tables.TableFormat; +import hep.dataforge.tables.ValueMap; +import hep.dataforge.values.Value; +import hep.dataforge.values.ValueType; import hep.dataforge.values.Values; import inr.numass.data.api.NumassBlock; import inr.numass.data.api.SignalProcessor; import org.jetbrains.annotations.Nullable; +import java.util.Map; + /** * An analyzer dispatcher which uses different analyzer for different meta * Created by darksnake on 11.07.2017. @@ -27,8 +34,7 @@ public class SmartAnalyzer extends AbstractAnalyzer { @Override public Values analyze(NumassBlock block, Meta config) { - //TODO add result caching - if(config.hasValue("type")) { + if (config.hasValue("type")) { switch (config.getString("type")) { case "simple": return simpleAnalyzer.analyze(block, config); @@ -40,11 +46,46 @@ public class SmartAnalyzer extends AbstractAnalyzer { throw new IllegalArgumentException("Analyzer not found"); } } else { - if(config.hasValue("t0")){ - return timeAnalyzer.analyze(block,config); + int t0 = getT0(block, config); + if (t0 == 0) { + Map map = simpleAnalyzer.analyze(block, config).asMap(); + map.putIfAbsent(TimeAnalyzer.T0_KEY, Value.of(0d)); + return new ValueMap(map); } else { - return simpleAnalyzer.analyze(block,config); + return timeAnalyzer.analyze(block, config.getBuilder().putValue(TimeAnalyzer.T0_KEY, t0)); } } } + + private double estimateCountRate(NumassBlock block) { + return (double) block.getEvents().count() / block.getLength().toMillis() * 1000; + } + + @ValueDef(name = "t0", type = ValueType.NUMBER, info = "Constant t0 cut") + @ValueDef(name = "t0.crFraction", type = ValueType.NUMBER, info = "The relative fraction of events that should be removed by time cut") + @ValueDef(name = "t0.min", type = ValueType.NUMBER, def = "0", info = "Minimal t0") + private int getT0(NumassBlock block, Meta meta) { + if (meta.hasValue("t0")) { + return meta.getInt("t0"); + } else if (meta.hasMeta("t0")) { + double fraction = meta.getDouble("t0.crFraction"); + double cr = estimateCountRate(block); + if (cr < meta.getDouble("t0.minCR", 0)) { + return 0; + } else { + return (int) Math.max(-1e9 / cr * Math.log(1d - fraction), meta.getDouble("t0.min", 0)); + } + } else { + return 0; + } + + } + + @Override + protected TableFormat getTableFormat(Meta config) { + if (config.hasValue(TimeAnalyzer.T0_KEY) || config.hasMeta(TimeAnalyzer.T0_KEY)) { + return timeAnalyzer.getTableFormat(config); + } + return super.getTableFormat(config); + } } diff --git a/numass-core/src/main/java/inr/numass/data/analyzers/TimeAnalyzer.java b/numass-core/src/main/java/inr/numass/data/analyzers/TimeAnalyzer.java index 77e11b3f..53b659c9 100644 --- a/numass-core/src/main/java/inr/numass/data/analyzers/TimeAnalyzer.java +++ b/numass-core/src/main/java/inr/numass/data/analyzers/TimeAnalyzer.java @@ -1,6 +1,8 @@ package inr.numass.data.analyzers; import hep.dataforge.meta.Meta; +import hep.dataforge.tables.TableFormat; +import hep.dataforge.tables.TableFormatBuilder; import hep.dataforge.tables.ValueMap; import hep.dataforge.values.Values; import inr.numass.data.api.NumassBlock; @@ -14,11 +16,18 @@ import java.util.concurrent.atomic.AtomicLong; import java.util.concurrent.atomic.AtomicReference; import java.util.stream.Stream; +import static hep.dataforge.tables.XYAdapter.*; +import static inr.numass.data.api.NumassPoint.HV_KEY; + /** * An analyzer which uses time information from events * Created by darksnake on 11.07.2017. */ public class TimeAnalyzer extends AbstractAnalyzer { + public static String T0_KEY = "t0"; + + public static String[] NAME_LIST = {LENGTH_KEY, COUNT_KEY, COUNT_RATE_KEY, COUNT_RATE_ERROR_KEY, WINDOW_KEY, TIME_KEY, T0_KEY}; + public static String[] NAME_LIST_WITH_HV = {HV_KEY, LENGTH_KEY, COUNT_KEY, COUNT_RATE_KEY, COUNT_RATE_ERROR_KEY, WINDOW_KEY, TIME_KEY, T0_KEY}; public TimeAnalyzer(@Nullable SignalProcessor processor) { super(processor); @@ -37,6 +46,7 @@ public class TimeAnalyzer extends AbstractAnalyzer { AtomicLong totalT = new AtomicLong(0); getEventsWithDelay(block, config) + .filter(pair -> pair.getValue() >= t0) .forEach(pair -> { totalN.incrementAndGet(); //TODO add progress listener here @@ -45,8 +55,9 @@ public class TimeAnalyzer extends AbstractAnalyzer { double countRate = 1e6 * totalN.get() / (totalT.get() / 1000 - t0 * totalN.get() / 1000);//1e9 / (totalT.get() / totalN.get() - t0); double countRateError = countRate / Math.sqrt(totalN.get()); - long count = (long) (totalT.get() * (countRate / 1e9)); - double length = totalT.get(); + double length = totalT.get() / 1e9; + long count = (long) (length * countRate); + if (block instanceof NumassPoint) { return ValueMap.of(NAME_LIST_WITH_HV, @@ -56,7 +67,9 @@ public class TimeAnalyzer extends AbstractAnalyzer { countRate, countRateError, new Integer[]{loChannel, upChannel}, - block.getStartTime()); + block.getStartTime(), + (double)t0 / 1000d + ); } else { return ValueMap.of(NAME_LIST, length, @@ -64,13 +77,14 @@ public class TimeAnalyzer extends AbstractAnalyzer { countRate, countRateError, new Integer[]{loChannel, upChannel}, - block.getStartTime() + block.getStartTime(), + (double)t0 / 1000d ); } } private long getT0(NumassBlock block, Meta config) { - return config.getValue("t0",0).longValue(); + return config.getValue("t0", 0).longValue(); } /** @@ -81,37 +95,46 @@ public class TimeAnalyzer extends AbstractAnalyzer { * @return */ public Stream> getEventsWithDelay(NumassBlock block, Meta config) { - long t0 = getT0(block, config); - AtomicReference lastEvent = new AtomicReference<>(null); Stream eventStream = super.getEvents(block, config);//using super implementation return eventStream.map(event -> { - if (lastEvent.get() == null) { - lastEvent.set(event); - return new Pair<>(event, 0L); - } else { - long res = event.getTimeOffset() - lastEvent.get().getTimeOffset(); - if (res >= 0) { - lastEvent.set(event); - return new Pair<>(event, res); - } else { - lastEvent.set(null); - return new Pair<>(event, 0L); - } + long res = lastEvent.get() == null ? -1L : event.getTimeOffset() - lastEvent.get().getTimeOffset(); + + if (res < 0) { + res = 0L; } - }).filter(pair -> pair.getValue() >= t0); + + lastEvent.set(event); + return new Pair<>(event, res); + }); } /** * The filtered stream of events + * * @param block * @param config * @return */ @Override public Stream getEvents(NumassBlock block, Meta config) { - return getEventsWithDelay(block, config).map(Pair::getKey); + long t0 = getT0(block, config); + return getEventsWithDelay(block, config).filter(pair -> pair.getValue() >= t0).map(Pair::getKey); + } + + @Override + protected TableFormat getTableFormat(Meta config) { + return new TableFormatBuilder() + .addNumber(HV_KEY, X_VALUE_KEY) + .addNumber(LENGTH_KEY) + .addNumber(COUNT_KEY) + .addNumber(COUNT_RATE_KEY, Y_VALUE_KEY) + .addNumber(COUNT_RATE_ERROR_KEY, Y_ERROR_KEY) + .addColumn(WINDOW_KEY) + .addTime() + .addNumber(T0_KEY) + .build(); } } diff --git a/numass-main/src/main/groovy/inr/numass/scripts/times/TestPointAnalyzer.groovy b/numass-main/src/main/groovy/inr/numass/scripts/times/TestPointAnalyzer.groovy index 46796b8f..907020bb 100644 --- a/numass-main/src/main/groovy/inr/numass/scripts/times/TestPointAnalyzer.groovy +++ b/numass-main/src/main/groovy/inr/numass/scripts/times/TestPointAnalyzer.groovy @@ -25,20 +25,20 @@ ctx.pluginManager().load(NumassPlugin.class) new GrindShell(ctx).eval { PlotHelper plot = plots - File rootDir = new File("D:\\Work\\Numass\\data\\2017_05\\Fill_3") + File rootDir = new File("D:\\Work\\Numass\\data\\2017_05\\Fill_2") NumassStorage storage = NumassStorageFactory.buildLocal(rootDir); - def set = "set_43" - def hv = 16000; + def set = "set_2" + def hv = 18300; def loader = storage.provide("loader::$set", NumassSet.class).get(); def point = loader.provide("$hv", NumassPoint.class).get() - def loChannel = 500; - def upChannel = 2000; + def loChannel = 450; + def upChannel = 3100; - def histogram = PointAnalyzer.histogram(point, loChannel, upChannel, 0.7, 1000).asTable(); + def histogram = PointAnalyzer.histogram(point, loChannel, upChannel, 1, 500).asTable(); println "finished histogram calculation..." diff --git a/numass-main/src/main/java/inr/numass/actions/MergeDataAction.java b/numass-main/src/main/java/inr/numass/actions/MergeDataAction.java index 1cfcb9f7..e82a3f64 100644 --- a/numass-main/src/main/java/inr/numass/actions/MergeDataAction.java +++ b/numass-main/src/main/java/inr/numass/actions/MergeDataAction.java @@ -99,16 +99,16 @@ public class MergeDataAction extends ManyToOneAction { return dp1; } - double Uset = dp1.getValue(parnames[0]).doubleValue(); + double voltage = dp1.getValue(NumassPoint.HV_KEY).doubleValue(); //усредняем измеренное напряжение - double Uread = (dp1.getValue(parnames[1]).doubleValue() + dp2.getValue(parnames[1]).doubleValue()) / 2; +// double Uread = (dp1.getValue(parnames[1]).doubleValue() + dp2.getValue(parnames[1]).doubleValue()) / 2; double t1 = dp1.getValue(NumassPoint.LENGTH_KEY).doubleValue(); double t2 = dp2.getValue(NumassPoint.LENGTH_KEY).doubleValue(); double time = t1 + t2; - long total = dp1.getValue(parnames[3]).intValue() + dp2.getValue(parnames[3]).intValue(); - long wind = dp1.getValue(parnames[4]).intValue() + dp2.getValue(parnames[4]).intValue(); + long total = dp1.getValue(NumassAnalyzer.COUNT_KEY).intValue() + dp2.getValue(NumassAnalyzer.COUNT_KEY).intValue(); +// long wind = dp1.getValue(parnames[4]).intValue() + dp2.getValue(parnames[4]).intValue(); // double corr = dp1.getValue(parnames[5]).doubleValue() + dp2.getValue(parnames[5]).doubleValue(); double cr1 = dp1.getValue(NumassAnalyzer.COUNT_RATE_KEY).doubleValue(); @@ -122,7 +122,7 @@ public class MergeDataAction extends ManyToOneAction { // абсолютные ошибки складываются квадратично double crErr = Math.sqrt(err1 * err1 * t1 * t1 + err2 * err2 * t2 * t2) / time; - ValueMap.Builder map = ValueMap.of(parnames, Uset, time, total, cr, crErr).builder(); + ValueMap.Builder map = ValueMap.of(parnames, voltage, time, total, cr, crErr).builder(); // if (dp1.getNames().contains("relCR") && dp2.getNames().contains("relCR")) { // double relCR = (dp1.getDouble("relCR") + dp2.getDouble("relCR")) / 2; @@ -141,7 +141,7 @@ public class MergeDataAction extends ManyToOneAction { throw new IllegalArgumentException(); } for (Values dp : d) { - double uset = dp.getValue(parnames[0]).doubleValue(); + double uset = dp.getValue(NumassPoint.HV_KEY).doubleValue(); if (!points.containsKey(uset)) { points.put(uset, new ArrayList<>()); } 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 aafab6d2..7459de04 100644 --- a/numass-main/src/main/java/inr/numass/models/LossCalculator.java +++ b/numass-main/src/main/java/inr/numass/models/LossCalculator.java @@ -53,6 +53,8 @@ public class LossCalculator { private final Map cache = new HashMap<>(); + + private LossCalculator() { cache.put(1, getSingleScatterFunction()); // cache.put(2, getDoubleScatterFunction()); @@ -386,7 +388,7 @@ public class LossCalculator { return integrator.integrate(integrand, 5d, margin); }; - return FunctionCaching.cacheUnivariateFunction(res, 0, margin, 200); + return FunctionCaching.cacheUnivariateFunction(0, margin, 200, res); } diff --git a/numass-main/src/main/java/inr/numass/models/sterile/NumassTransmission.java b/numass-main/src/main/java/inr/numass/models/sterile/NumassTransmission.java index fc6319c1..d0992160 100644 --- a/numass-main/src/main/java/inr/numass/models/sterile/NumassTransmission.java +++ b/numass-main/src/main/java/inr/numass/models/sterile/NumassTransmission.java @@ -13,13 +13,13 @@ import hep.dataforge.values.Values; import inr.numass.models.LossCalculator; import inr.numass.utils.ExpressionUtils; import org.apache.commons.math3.analysis.BivariateFunction; +import org.apache.commons.math3.analysis.UnivariateFunction; import org.slf4j.LoggerFactory; import java.util.HashMap; import java.util.Map; /** - * * @author Alexander Nozik */ public class NumassTransmission extends AbstractParametricBiFunction { @@ -28,13 +28,14 @@ public class NumassTransmission extends AbstractParametricBiFunction { private static final double ION_POTENTIAL = 15.4;//eV private final LossCalculator calculator; private final BivariateFunction trapFunc; + private Map lossCache = new HashMap<>(); private final boolean adjustX; public NumassTransmission(Context context, Meta meta) { super(list); this.calculator = LossCalculator.instance(); - adjustX = meta.getBoolean("adjustX",false); + adjustX = meta.getBoolean("adjustX", false); if (meta.hasValue("trapping")) { String trapFuncStr = meta.getString("trapping"); if (trapFuncStr.startsWith("function::")) { @@ -54,7 +55,7 @@ public class NumassTransmission extends AbstractParametricBiFunction { } public double getX(double eIn, Values set) { - if(adjustX){ + if (adjustX) { //From our article return set.getDouble("X") * Math.log(eIn / ION_POTENTIAL) * eIn * ION_POTENTIAL / 1.9580741410115568e6; } else { @@ -93,6 +94,18 @@ public class NumassTransmission extends AbstractParametricBiFunction { double X = getX(eIn, set); // loss part double loss = calculator.getTotalLossValue(X, eIn, eOut); +// double loss; +// +// if(eIn-eOut >= 300){ +// loss = 0; +// } else { +// UnivariateFunction lossFunction = this.lossCache.computeIfAbsent(X, theX -> +// FunctionCaching.cacheUnivariateFunction(0, 300, 400, x -> calculator.getTotalLossValue(theX, eIn, eIn - x)) +// ); +// +// loss = lossFunction.value(eIn - eOut); +// } + //trapping part double trap = getParameter("trap", set) * trapFunc.value(eIn, eOut); return loss + trap;