diff --git a/eos/event_data_types.py b/eos/event_data_types.py index d480ba8..466d463 100644 --- a/eos/event_data_types.py +++ b/eos/event_data_types.py @@ -39,6 +39,7 @@ EVENT_TYPE = np.dtype([('tof', np.float64), ('pixelID', np.uint32), ('mask', np. PACKET_TYPE = np.dtype([('start_index', np.uint32), ('time', np.int64)]) PULSE_TYPE = np.dtype([('time', np.int64), ('monitor', np.float32)]) PC_TYPE = np.dtype([('current', np.float32), ('time', np.int64)]) +BM_TYPE = np.dtype([('value', np.int64), ('time', np.int64)]) LOG_TYPE = np.dtype([('value', np.float32), ('time', np.int64)]) # define the bitmask for individual event filters @@ -68,6 +69,7 @@ class AmorEventStream: packets: np.recarray # PACKET_TYPE pulses: np.recarray # PULSE_TYPE proton_current: np.recarray # PC_TYPE + beam_monitor: np.recarray # BM_TYPE device_logs: Dict[str, np.recarray] = field(default_factory=dict) # LOG_TYPE class EventDatasetProtocol(Protocol): diff --git a/eos/event_handling.py b/eos/event_handling.py index 84c70c3..d062239 100644 --- a/eos/event_handling.py +++ b/eos/event_handling.py @@ -91,6 +91,8 @@ class AssociatePulseWithMonitor(EventDataAction): * 2*dataset.timing.tau * 1e-3 elif self.monitorType==MonitorType.time: dataset.data.pulses.monitor = 2*dataset.timing.tau + elif self.monitorType==MonitorType.neutron_monitor: + dataset.data.pulses.monitor = dataset.data.beam_monitor.value / 1e4 else: # pulses dataset.data.pulses.monitor = 1 @@ -201,4 +203,4 @@ class ApplyMask(EventDataAction): d.device_logs[key] = np.recarray(d.pulses.shape, dtype = log.dtype) d.device_logs[key].time = d.pulses.time d.device_logs[key].value = d.pulses[key] - dataset.update_info_from_logs() \ No newline at end of file + dataset.update_info_from_logs() diff --git a/eos/file_reader.py b/eos/file_reader.py index ad30f68..f6d1fcc 100644 --- a/eos/file_reader.py +++ b/eos/file_reader.py @@ -7,6 +7,7 @@ import sys import h5py import numpy as np import logging +import subprocess from datetime import datetime @@ -16,7 +17,7 @@ from orsopy.fileio.model_language import SampleModel from . import const, compat from .header import Header from .event_data_types import AmorGeometry, AmorTiming, AmorEventStream, LOG_TYPE, PACKET_TYPE, EVENT_TYPE, PULSE_TYPE, \ - PC_TYPE, DA00EventBuffer + PC_TYPE, BM_TYPE, DA00EventBuffer try: import zoneinfo @@ -74,7 +75,7 @@ class AmorHeader: data=('entry1/detector/data', None), # data group used to load events from trigger=('entry1/detector/trigger', float), - monitor=('entry1/beam_monitor', None), # special data group (da00 flat-buffer scheme) + beam_monitor=('entry1/beam_monitor', None), # special data group (da00 flat-buffer scheme) proton_current=('entry1/detector/proton_current', float), acquisition_filter=('entry1/detector/acquisition_filter', int), ) @@ -457,14 +458,15 @@ class AmorEventData(AmorHeader): events.pixelID = data_group['event_id'][self.first_index:self.last_index+1] events.mask = 0 - self._monitor = self.rv('monitor', default=None) + self._beam_monitor = self.rv('beam_monitor', default=None) pulses = self.read_chopper_trigger_stream(packets) current = self.read_proton_current_stream(packets) + beam_monitor = self.read_proton_current_stream(packets) # read parameter logs not present in old files to ensure they are in self._log_keys if they exist _acquisition_filter = self.rv('acquisition_filter', default=None) - self.data = AmorEventStream(events, packets, pulses, current) + self.data = AmorEventStream(events, packets, pulses, current, beam_monitor) if self.first_index>0 and not self.EOF: # label the file name if not all events were used @@ -521,8 +523,8 @@ class AmorEventData(AmorHeader): pulseTimeS = chopper1TriggerTime else: logging.critical(' No chopper trigger data available, using event steram instead, pulse filtering will fail!') - startTime = np.array(self._monitor.time[0], dtype=np.int64) - pulseTimeS = self._monitor.time + startTime = np.array(self._beam_monitor.time[0], dtype=np.int64) + pulseTimeS = self._beam_monitor.time pulses = np.recarray(pulseTimeS.shape, dtype=PULSE_TYPE) pulses.time = pulseTimeS pulses.monitor = 1. # default is monitor pulses as it requires no calculation @@ -546,6 +548,20 @@ class AmorEventData(AmorHeader): (proton_current.time<=packets.time[-1])] return proton_current + def read_beam_monitor_stream(self, packets): + _signal = self._beam_monitor.signal + _time = self._beam_monitor.time + plse = np.array([0, 5, 6, 11]) + base = np.array([2, 3, 8, 9]) + i = 0 + beam_monitor = np.recarray(shape=np.shape(packets.time)[0], dtype=BM_TYPE) + for tt, tme in enumerate(_time): + if (tme>=packets.time[0])&(tme<=packets.time[-1]): + beam_monitor.value[i] = _signal[tt][plse].sum() - _signal[tt][base].sum() + beam_monitor.time[i] = tme + i += 1 + return beam_monitor[:i+1] + def info(self): output = "" for key in ['owner', 'experiment', 'sample', 'instrument_settings']: @@ -562,10 +578,11 @@ class AmorEventData(AmorHeader): new_events = np.concatenate([self.data.events, other.data.events]).view(np.recarray) new_pulses = np.concatenate([self.data.pulses, other.data.pulses]).view(np.recarray) new_proton_current = np.concatenate([self.data.proton_current, other.data.proton_current]).view(np.recarray) + new_beam_monitor = np.concatenate([self.data.beam_monitor, other.data.beam_monitor]).view(np.recarray) new_packets = np.concatenate([self.data.packets, other.data.packets]).view(np.recarray) new_packets.start_index[self.data.packets.shape[0]:] += self.data.events.shape[0] - self.data = AmorEventStream(new_events, new_packets, new_pulses, new_proton_current) - # Indicate that this is amodified dataset, basically counts number of appends as negative indices + self.data = AmorEventStream(new_events, new_packets, new_pulses, new_proton_current, new_beam_monitor) + # Indicate that this is a modified dataset, basically counts number of appends as negative indices self.last_index = min(self.last_index-1, -1) self.file_list += other.file_list @@ -593,5 +610,6 @@ class AmorEventData(AmorHeader): setattr(output, key, value) # TODO: this is not strictly correct, as the packet/event relationship is lost output.data = AmorEventStream(self.data.events[event_filter], self.data.packets, - self.data.pulses[pulse_filter], self.data.proton_current) + self.data.pulses[pulse_filter], self.data.proton_current, + self.data.beam_monitor) return output diff --git a/eos/options.py b/eos/options.py index 4eabd24..10db230 100644 --- a/eos/options.py +++ b/eos/options.py @@ -247,7 +247,7 @@ class MonitorType(StrEnum): debug = 'x' MONITOR_UNITS = { - MonitorType.neutron_monitor: 'cnts', + MonitorType.neutron_monitor: 'au', MonitorType.proton_charge: 'mC', MonitorType.time: 's', MonitorType.auto: 'various', @@ -273,7 +273,7 @@ class ExperimentConfig(ArgParsable): }, ) chopperSpeed: float = field( - default=500, + default=450, metadata={ 'short': 'cs', 'group': 'instrument settings', @@ -434,7 +434,7 @@ class ReflectivityReductionConfig(ArgParsable): 'group': 'input data', 'help': 'apply convolution on lambda axes to smooth the normalization data, useful for low statistics'}) scale: List[float] = field( - default_factory=lambda: [1.], + default_factory=lambda: [], metadata={ 'short': 's', 'group': 'data manicure', diff --git a/eos/projection.py b/eos/projection.py index 2cdc50e..c9db818 100644 --- a/eos/projection.py +++ b/eos/projection.py @@ -248,9 +248,9 @@ class LZProjection(ProjectionInterface): def scale(self, factor: float): if not self.is_normalized: raise ValueError("Dataset needs to be normalized, first") - self.data.ref *= factor - self.data.err *= factor - self.norm_monitor /= factor + self.data.ref /= factor + self.data.err /= factor + self.norm_monitor *= factor def project_on_qz(self): if not self.is_normalized: @@ -829,4 +829,4 @@ class CombinedProjection(ProjectionInterface): pi.update_plot() # def update_range(self, event): - # ... \ No newline at end of file + # ... diff --git a/eos/reduction_reflectivity.py b/eos/reduction_reflectivity.py index ce8fc3b..6168c38 100644 --- a/eos/reduction_reflectivity.py +++ b/eos/reduction_reflectivity.py @@ -204,11 +204,12 @@ class ReflectivityReduction: logging.info(f' monitor = {self.monitor:8.2f} {MONITOR_UNITS[self.config.experiment.monitorType]}') proj:LZProjection = self.project_on_lz() - try: - scale = self.config.reduction.scale[i] - except IndexError: - scale = self.config.reduction.scale[-1] - proj.scale(scale) + if not self.config.reduction.is_default('scale'): + try: + scale = self.config.reduction.scale[i] + except IndexError: + scale = self.config.reduction.scale[-1] + proj.scale(scale) if 'Rqz.ort' in self.config.output.outputFormats: headerRqz = self.header.orso_header()