From 44ab8fb211534454e02389b57c12fbc1405c7c3f Mon Sep 17 00:00:00 2001 From: Artur Glavic Date: Mon, 4 Mar 2024 15:26:29 +0100 Subject: [PATCH] Start refactoring reduction into multiple methods --- libeos/dataset.py | 10 +- libeos/reduction.py | 380 ++++++++++++++++++++++---------------------- 2 files changed, 201 insertions(+), 189 deletions(-) diff --git a/libeos/dataset.py b/libeos/dataset.py index 6c06fe4..117f509 100644 --- a/libeos/dataset.py +++ b/libeos/dataset.py @@ -13,7 +13,6 @@ from . import __version__, const from .instrument import Detector from .options import ExperimentConfig, ReaderConfig - class Header: """orso compatible output file header content""" @@ -58,6 +57,15 @@ class Header: ] return cols + def orso_header(self, columns=None, extra_columns=[]): + """ + Generate ORSO header from a copy of this class' data. + """ + ds = fileio.DataSource.from_dict(self.data_source().to_dict()) + red = fileio.Reduction.from_dict(self.reduction.to_dict()) + if columns is None: + columns = self.columns() + return fileio.Orso(ds, red, columns+extra_columns) class AmorData: """read meta-data and event streams from .hdf file(s), apply filters and conversions""" diff --git a/libeos/reduction.py b/libeos/reduction.py index 1f62511..938e78c 100644 --- a/libeos/reduction.py +++ b/libeos/reduction.py @@ -24,197 +24,202 @@ class AmorReduction: if not os.path.exists(f'{self.reader_config.dataPath}'): logging.debug(f'Creating destination path {self.reader_config.dataPath}') os.system(f'mkdir {self.reader_config.dataPath}') - fromHDF = AmorData(self.startTime, header=self.header, reader_config=self.reader_config, config=self.experiment_config) # load or create normalisation matrix if self.reduction_config.normalisationFileIdentifier: - normalise = True - norm_lz, normAngle, normFileList = self.normalisation_map(self.reduction_config.normalisationFileIdentifier[0]) - self.header.reduction.corrections.append('normalisation with \'additional files\'') + self.create_normalisation_map(self.reduction_config.normalisationFileIdentifier[0]) else: - normalise = False - norm_lz = self.grid.lz() - normAngle = 1. + self.norm_lz = self.grid.lz() + self.normAngle = 1. logging.warning('normalisation matrix: none requested') # load R(q_z) curve to be subtracted: if self.reduction_config.subtract: - sq_q, sR_q, sdR_q, sFileName = self.loadRqz(self.reduction_config.subtract) - subtract = True - logging.warning(f'loaded background file: {sFileName}') - self.header.reduction.corrections.append(f'background from \'{sFileName}\' subtracted') + self.sq_q, self.sR_q, self.sdR_q, self.sFileName = self.loadRqz(self.reduction_config.subtract) + logging.warning(f'loaded background file: {self.sFileName}') + self.header.reduction.corrections.append(f'background from \'{self.sFileName}\' subtracted') + self.subtract = True else: - subtract = False + self.subtract = False # load measurement data and do the reduction - datasetsRqz = [] - datasetsRlt = [] + self.datasetsRqz = [] + self.datasetsRlt = [] + self.file_reader = AmorData(self.startTime, header=self.header, reader_config=self.reader_config, config=self.experiment_config) for i, short_notation in enumerate(self.reduction_config.fileIdentifier): - logging.warning('reading input:') - self.header.measurement_data_files = [] - fromHDF.read_data(short_notation) - - if self.reduction_config.timeSlize: - wallTime_e = fromHDF.wallTime_e - columns = self.header.columns()+[fileio.Column('time', 's', 'time relative to start of measurement series')] - headerRqz = fileio.Orso(self.header.data_source(), self.header.reduction, columns) - - interval = self.reduction_config.timeSlize[0] - try: - start = self.reduction_config.timeSlize[1] - except: - start = 0 - try: - stop = self.reduction_config.timeSlize[2] - except: - stop = wallTime_e[-1] - # make overwriting log lines possible by removing newline at the end - logging.StreamHandler.terminator = "\r" - for i, time in enumerate(np.arange(start, stop, interval)): - logging.warning(f' time slize {i:4d}') - - filter_e = np.where((time0, norm_lz, np.nan) + self.norm_lz, bins_l, bins_z = np.histogram2d(lamda_e, detZ_e, bins = (self.grid.lamda(), self.grid.z())) + self.norm_lz = np.where(self.norm_lz>0, self.norm_lz, np.nan) # correct for the SM reflectivity lamda_l = self.grid.lamda() - theta_z = normAngle + fromHDF.delta_z + theta_z = self.normAngle + fromHDF.delta_z lamda_lz = (self.grid.lz().T*lamda_l[:-1]).T theta_lz = self.grid.lz()*theta_z qz_lz = 4.0*np.pi * np.sin(np.deg2rad(theta_lz)) / lamda_lz Rsm_lz = np.ones(np.shape(qz_lz)) Rsm_lz = np.where(qz_lz>0.0217, 1-(qz_lz-0.0217)*(0.0625/0.0217), Rsm_lz) Rsm_lz = np.where(qz_lz>0.0217*5, np.nan, Rsm_lz) - norm_lz = norm_lz / Rsm_lz + self.norm_lz = self.norm_lz / Rsm_lz if len(lamda_e) > 1e6: head = ('normalisation matrix based on the measurements\n' f'{fromHDF.file_list}\n' - f'nu - mu = {normAngle}\n' - f'shape= {np.shape(norm_lz)} (lambda, z)\n' + f'nu - mu = {self.normAngle}\n' + f'shape= {np.shape(self.norm_lz)} (lambda, z)\n' f'measured at mu = {fromHDF.mu:6.3f} deg\n' f'N(l_lambda, z) = theta(z) / sum_i=-1..1 I(l_lambda+i, z)') head = head.replace('../', '') head = head.replace('./', '') head = head.replace('raw/', '') - np.savetxt(f'{dataPath}/{name}.norm', norm_lz, header = head) - normFileList = fromHDF.file_list - return norm_lz, normAngle, normFileList - + np.savetxt(f'{dataPath}/{name}.norm', self.norm_lz, header = head) + self.normFileList = fromHDF.file_list + self.header.reduction.corrections.append('normalisation with \'additional files\'') def project_on_lz(self, fromHDF, norm_lz, normAngle, lamda_e, detZ_e): # projection on lambda-z-grid