diff --git a/src/backscattering_analyzer/__init__.py b/src/backscattering_analyzer/__init__.py index 2a72377..2ee35dc 100644 --- a/src/backscattering_analyzer/__init__.py +++ b/src/backscattering_analyzer/__init__.py @@ -7,12 +7,16 @@ install(__name__) logger = getLogger(__name__) -def show_projection(experiment: "Experiment") -> None: +def show_projection( + experiment: "Experiment", + start_frequency: float | None = None, + end_frequency: float | None = None, +) -> None: """ Show projection data with matplotlib """ logger.debug(_("showing experiment result")) - from matplotlib.pyplot import loglog, show, legend, close + from matplotlib.pyplot import loglog, show, legend, close, vlines excited = experiment.signals["excited"].psd().sqrt() reference = experiment.signals["reference"].psd().sqrt() @@ -28,6 +32,24 @@ def show_projection(experiment: "Experiment") -> None: (experiment.projection + reference).y, label="sum reference + excited", ) # type: ignore[ReportUnusedCallResult] + if start_frequency is not None: + vlines( + [ + start_frequency, + ], + min(experiment.projection.y), + max(experiment.projection.y), + color = 'k', + ) # type: ignore[ReportUnusedCallResult] + if end_frequency is not None: + vlines( + [ + end_frequency, + ], + min(experiment.projection.y), + max(experiment.projection.y), + color = 'k', + ) # type: ignore[ReportUnusedCallResult] legend() # type: ignore[ReportUnusedCallResult] show() close() diff --git a/src/backscattering_analyzer/experiment.py b/src/backscattering_analyzer/experiment.py index 799b9f9..f6a7661 100644 --- a/src/backscattering_analyzer/experiment.py +++ b/src/backscattering_analyzer/experiment.py @@ -1,7 +1,7 @@ from pathlib import Path from tomllib import load -from numpy import argmin, float64, pi +from numpy import argmin, float64, pi, log10 from numpy.core.function_base import logspace from numpy.core.multiarray import array from numpy.typing import NDArray @@ -295,9 +295,11 @@ class Experiment: def fit_factors( self, - start: float = 15, - end: float = 100, - x: NDArray[float64] | None = None, + scatter_min: float = 1e-20, + scatter_max: float = 1, + start_frequency: float = 15, + end_frequency: float = 100, + precision: int = 3, ) -> dict[str, float]: """ Find the best factor (first order only) to get the projection @@ -311,35 +313,51 @@ class Experiment: excited, reference, phase, - ) = self.get_factors(start=start, end=end) + ) = self.get_factors(start=start_frequency, end=end_frequency) - if x is None: - x = logspace(-10, 0, 1000) + for index in range(precision): + logger.debug( + _("search for a local minimum") + ) - 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)) - ] - ) + x = logspace(log10(scatter_min), log10(scatter_max), 1000) - pre_scatter_factor: float = x[argmin(y)] + 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)) + ] + ) + + if argmin(y) == 0: + logger.warning(_("smaller than current range allows")) + scatter_max: float = x[1] + elif argmin(y) == len(x) - 1: + logger.warning(_("bigger than current range allows")) + scatter_min: float = x[-2] + else: + scatter_min: float = x[argmin(y) - 1] + scatter_max: float = x[argmin(y) + 1] + + logger.debug(_("local minimum found")) + + pre_scatter_factor: float = scatter_min / 2 + scatter_max / 2 logger.info( _("found a scattering factor of {factor}").format(