From 2bb202a622bbc1fb56f3a023394c649bab538699 Mon Sep 17 00:00:00 2001 From: linarphy Date: Thu, 3 Apr 2025 13:33:53 +0200 Subject: [PATCH] Remove old code --- src/backscattering_analyzer/experiment.py.old | 461 ------------------ 1 file changed, 461 deletions(-) delete mode 100644 src/backscattering_analyzer/experiment.py.old diff --git a/src/backscattering_analyzer/experiment.py.old b/src/backscattering_analyzer/experiment.py.old deleted file mode 100644 index 4678678..0000000 --- a/src/backscattering_analyzer/experiment.py.old +++ /dev/null @@ -1,461 +0,0 @@ -from pathlib import Path -from tomllib import load - -from numpy import argmin, pi, log10, logspace, array -from science_signal import interpolate -from science_signal.signal import Signal -from backscattering_analyzer import logger -from backscattering_analyzer.acquisition import AcquisitionType -from backscattering_analyzer.data import ( - Data, - load_coupling, -) -from backscattering_analyzer.utils import ( - compute_light, - process_movement, -) - - -class Experiment: - """ - An experiment in Virgo consisting of exciting (oscillating) a bench - """ - - def __init__( - self, - bench: str, - date: str, - values_file: str, - vibration_frequency: float = 0.2, - ) -> None: - logger.info( - _("initializing the experiment on the bench {bench} in {date}").format( - bench=bench, date=date - ) - ) - self.bench = bench - self.date = date - self.vibration_frequency = vibration_frequency - try: - with open(values_file, "rb") as file: - values = load(file) - logger.debug( - _("successfully loaded values from {file}").format(file=values_file) - ) - except OSError: - logger.critical(_("{file} cannot be read").format(file=values_file)) - raise Exception(_("{file] cannot be read").format(file=values_file)) - - try: - self.wavelength: float = values["interferometer"]["wavelength"] - except KeyError: - logger.warning( - _( - "cannot find the wavelength in {file}, default" - + "value will be used" - ).format(file=values_file) - ) - self.wavelength = 1.064e-6 - logger.debug( - _("value of wavelength: {wavelength}").format(wavelength=self.wavelength) - ) - self.calibrations: dict[str, float] = { - "bench": 1.0, - "mirror": 1.0, - } - try: - self.calibrations = { - "bench": values["benches"][bench]["calib"]["bench"], - "mirror": values["benches"][bench]["calib"]["mirror"], - } - except KeyError as e: - logger.error( - _( - "cannot find the calibration in {file}, default " - + "value will be used, the result will likely be " - + "wrong: {error}" - ).format(file=values_file, error=e) - ) - logger.debug( - _("value of calibrations (bench, mirror): ({bench}, " + "{mirror})").format( - bench=self.calibrations["bench"], - mirror=self.calibrations["mirror"], - ) - ) - self.factors: dict[str, float | None] = { - "pre": None, - "true": None, - } - try: - self.factors = { - "pre": 1e6 * values["benches"][bench]["scatter_factor"][0], - "true": values["benches"][bench]["scatter_factor"][0], - } - except KeyError: - logger.warning( - _( - "cannot find the scatter factor in {file}, no " - + "value will be set and the scatter factor will " - + "be searched with a fit" - ).format(file=values_file) - ) - logger.debug( - _( - "values of scattering factor (outscale, real): ({pre}, " + "{true})" - ).format( - pre=self.factors["pre"], - true=self.factors["true"], - ) - ) - self.powers: dict[str, float] = { - "laser": 23, - "dark_fringe": 8e-3, - } - try: - self.powers = { - "laser": values["dates"][date]["laser"], - "dark_fringe": values["dates"][date]["dark_fringe"], - } - except KeyError as e: - logger.error( - _( - "cannot find the power values in {file}, default " - + "one will be used, the result will likely be " - + "wrong: {error}" - ).format(file=values_file, error=e) - ) - logger.debug( - _( - "values of powers (laser, dark_fringe): ({laser}, " + "{dark_fringe})" - ).format( - laser=self.powers["laser"], - dark_fringe=self.powers["dark_fringe"], - ) - ) - try: - self.modelisation_file: str = values["dates"][date]["modelisation"] - except KeyError: - logger.error( - _( - "cannot find the modelisation file linked to this " - + "date in {file}, default one will be used, the " - + "result will likely be wrong" - ).format(file=values_file) - ) - self.modelisation_file = "scatterCouplingO4.mat" - logger.debug( - _("modelisation file: {file}").format( - file=self.modelisation_file, - ) - ) - data_folder_name: str = "/home/demagny/data" - try: - data_folder_name = values["general"]["data_folder"] - except KeyError: - logger.error( - _( - "cannot find the data folder containing the " - + "signals values in {file}, default one will be " - + "set, if its doesn't exist the program will crash" - ).format(file=values_file) - ) - logger.debug( - _("data folder containing signal values: {folder}").format( - folder=data_folder_name, - ) - ) - - data_folder = Path(data_folder_name) - - self.data = Data(data_folder, bench, date) - - self.reference_movement: None | Signal = process_movement( - self.data.reference.movements.bench, - self.data.reference.movements.mirror, - self.calibrations, - vibration_frequency, - ) - self.excited_movement: None | Signal = None - try: - self.excited_movement = process_movement( - self.data.excited.movements.bench, - self.data.excited.movements.mirror, - self.calibrations, - vibration_frequency, - ) - except AttributeError: - pass - logger.debug(_("movement processed")) - - model_powers: dict[str, float] = { - "laser": 0, - "dark_fringe": 0, - } - correct_power: bool = False - try: - correct_power = values["dates"][date]["correct-power"] - if correct_power: - model_powers = { - "laser": values["dates"][date]["coupling"]["laser"], - "dark_fringe": values["dates"][date]["coupling"]["dark_fringe"], - } - except KeyError as error: - print(error) - logger.warning( - _( - "cannot find if coupling values should be " - + "corrected due to the diff between power used in " - + "model and the real one. No correction will be " - + "applied" - ) - ) - - try: - self.coupling = load_coupling( - data_folder, - bench, - date, - self.modelisation_file, - self.powers, - model_powers, - ) - logger.debug(_("coupling signals successfully loaded")) - except OSError: - logger.critical(_("cannot load coupling data")) - raise Exception(_("cannot load coupling data")) - - if self.factors["pre"] is None: - logger.debug( - _( - "scattering factor is not defined, the script will " - + "try to find the right one by fitting" - ) - ) - self.factors = self.fit_factors() - - self.projection_reference = self.compute_projection( - AcquisitionType.REFERENCE, - ) - try: - self.projection_excited = self.compute_projection( - AcquisitionType.EXCITED, - ) - except Exception: - pass - logger.debug(_("end of experiment initialisation")) - - def get_factors( - self, - phase: float, - *signals: Signal, - start: float | None = None, - end: float | None = None, - type: AcquisitionType = AcquisitionType.EXCITED, - ) -> tuple[Signal, ...]: - """ - Get factor to compute the projection for a given scatter factor - """ - coupling_n = self.coupling[0] - coupling_d = self.coupling[1] - - if type is AcquisitionType.EXCITED: - if self.excited_movement is not None: - factor_n = (self.excited_movement * phase).sin().psd().sqrt() - factor_d = (self.excited_movement * phase).cos().psd().sqrt() - else: - logger.critical( - _( - "no excited signal given, cannot get factors " - + "of the excited signal" - ) - ) - raise Exception( - _( - "no excited signal given, cannot get factors " - + "of the excited signal" - ) - ) - elif type is AcquisitionType.REFERENCE: - if self.reference_movement is not None: - factor_n = (self.reference_movement * phase).sin().psd().sqrt() - factor_d = (self.reference_movement * phase).cos().psd().sqrt() - else: - logger.critical( - _( - "no reference signal given, cannot get factors " - + "of the reference signal" - ) - ) - raise Exception( - _( - "no reference signal given, cannot get factors " - + "of the reference signal" - ) - ) - else: - logger.critical(_("unknown acquisition type")) - raise Exception(_("unknown acquisition type")) - - if start is None: - start = coupling_n.x[0] - if end is None: - end = coupling_n.x[-1] - - coupling_n = coupling_n.cut(start, end) - - return interpolate( - factor_n, - coupling_n, - factor_d, - coupling_d, - *signals, - ) - - def fit_factors( - self, - scatter_min: float = 1e-20, - scatter_max: float = 1, - start_frequency: float = 20, - end_frequency: float = 35, - precision: int = 3, - ) -> dict[str, None | float]: - """ - Find the best factor (first order only) to get the projection - that best match the excited signal - """ - phase = 4 * pi / self.wavelength - ( - factor_n, - coupling_n, - factor_d, - coupling_d, - excited, - reference, - ) = self.get_factors( - phase, - self.data.excited.sensibility.psd().sqrt(), - self.data.reference.sensibility.psd().sqrt(), - start=start_frequency, - end=end_frequency, - type=AcquisitionType.EXCITED, - ) - - scatter_max = 1e10 - scatter_min = 1e-10 - - for index in range(precision): - logger.debug(_("search for a local minimum")) - - x = logspace(log10(scatter_min), log10(scatter_max), 1000) - - y = array( - [ - sum( - abs( - Signal( - factor_n.x, - compute_light( - pre_scatter_factor=x[i], - factor_n=factor_n, - coupling_n=coupling_n, - factor_d=factor_d, - coupling_d=coupling_d, - phase=phase, - ), - ) - + reference - - excited - ).y - ) - for i in range(len(x)) - ] - ) - from matplotlib.pyplot import figure, show - - diff_1 = abs( - Signal( - factor_n.x, - compute_light( - pre_scatter_factor=x[100], - factor_n=factor_n, - coupling_n=coupling_n, - factor_d=factor_d, - coupling_d=coupling_d, - phase=phase, - ), - ) - + reference - - excited - ) - diff_2 = abs( - Signal( - factor_n.x, - compute_light( - pre_scatter_factor=x[900], - factor_n=factor_n, - coupling_n=coupling_n, - factor_d=factor_d, - coupling_d=coupling_d, - phase=phase, - ), - ) - + reference - - excited - ) - - fig = figure() - fig.gca().loglog(diff_1.x, diff_1.y, label=str(x[100])) - fig.gca().loglog(diff_2.x, diff_2.y, label=str(x[900])) - fig.gca().legend() - show() - - if argmin(y) == 0: - logger.warning(_("smaller than current range allows")) - scatter_max = x[1] - elif argmin(y) == len(x) - 1: - logger.warning(_("bigger than current range allows")) - scatter_min = x[-2] - else: - scatter_min = x[argmin(y) - 1] - scatter_max = x[argmin(y) + 1] - - logger.debug(_("local minimum found")) - - pre_scatter_factor = scatter_min / 2 + scatter_max / 2 - - logger.info( - _("found a scattering factor of {factor}").format( - factor=pre_scatter_factor * 1e-6 - ) - ) - - return { - "pre": pre_scatter_factor, - "true": pre_scatter_factor * 1e-6, - } - - def compute_projection( - self, type: AcquisitionType = AcquisitionType.EXCITED - ) -> Signal: - """ - Compute the projection of the noise from current values - """ - phase = 4 * pi / self.wavelength - ( - factor_n, - coupling_n, - factor_d, - coupling_d, - ) = self.get_factors(phase, type=type) - - return Signal( - factor_n.x, - compute_light( - pre_scatter_factor=self.factors["pre"], - factor_n=factor_n, - coupling_n=coupling_n, - factor_d=factor_d, - coupling_d=coupling_d, - phase=phase, - ), - )