diff --git a/qudi_hira_analysis/_fitmethods/antibunchingmethods.py b/qudi_hira_analysis/_fitmethods/antibunchingmethods.py index 2777d53..a3481b2 100644 --- a/qudi_hira_analysis/_fitmethods/antibunchingmethods.py +++ b/qudi_hira_analysis/_fitmethods/antibunchingmethods.py @@ -5,7 +5,8 @@ def make_antibunching_model(self, prefix=None): - def antibunching(x: np.ndarray, n: float, a: float, b: float, tau0: float, tau1: float, tau2: float) \ + def antibunching(x: np.ndarray, n: float, a: float, b: float, tau0: float, + tau1: float, tau2: float) \ -> np.ndarray: """ Fit to function @@ -34,7 +35,8 @@ def estimate_antibunching_dip(self, x_axis, data, params): return error, params -def make_antibunching_fit(self, x_axis, data, estimator, units=None, add_params=None, **kwargs): +def make_antibunching_fit(self, x_axis, data, estimator, units=None, add_params=None, + **kwargs): model, params = self.make_antibunching_model() error, params = estimator(x_axis, data, params) @@ -47,7 +49,7 @@ def make_antibunching_fit(self, x_axis, data, estimator, units=None, add_params= except: result = model.fit(data, x=x_axis, params=params, **kwargs) self.log.warning('The 1D antibunching fit did not work. Error ' - 'message: {0}\n'.format(result.message)) + f'message: {result.message}\n') # Write the parameters to allow human-readable output to be generated result_str_dict = OrderedDict() diff --git a/qudi_hira_analysis/_fitmethods/decaylikemethods.py b/qudi_hira_analysis/_fitmethods/decaylikemethods.py index 0c66f01..84ebb6f 100644 --- a/qudi_hira_analysis/_fitmethods/decaylikemethods.py +++ b/qudi_hira_analysis/_fitmethods/decaylikemethods.py @@ -1,4 +1,3 @@ -# -*- coding: utf-8 -*- """ This file contains methods for decay-like fitting, these methods are imported by class FitLogic. @@ -71,7 +70,7 @@ def barestretchedexponentialdecay_function(x, beta, lifetime): if not isinstance(prefix, str) and prefix is not None: - self.log.error('The passed prefix <{0}> of type {1} is not a string and' + self.log.error('The passed prefix <{}> of type {} is not a string and' 'cannot be used as a prefix and will be ignored for now.' 'Correct that!'.format(prefix, type(prefix))) model = Model(barestretchedexponentialdecay_function, @@ -101,7 +100,8 @@ def make_bareexponentialdecay_model(self, prefix=None): the method make_barestretchedexponentialdecay_model. """ - bare_exp_decay, params = self.make_barestretchedexponentialdecay_model(prefix=prefix) + bare_exp_decay, params = self.make_barestretchedexponentialdecay_model( + prefix=prefix) bare_exp_decay.set_param_hint(name='beta', value=1, vary=False) params = bare_exp_decay.make_params() @@ -149,7 +149,8 @@ def make_decayexponentialstretched_model(self, prefix=None): the method make_barestretchedexponentialdecay_model. """ - bare_stre_exp_decay, params = self.make_barestretchedexponentialdecay_model(prefix=prefix) + bare_stre_exp_decay, params = self.make_barestretchedexponentialdecay_model( + prefix=prefix) amplitude_model, params = self.make_amplitude_model() constant_model, params = self.make_constant_model(prefix=prefix) @@ -199,7 +200,8 @@ def make_biexponential_model(self, prefix=None): # single exponential decay with offset # ########################################## -def make_decayexponential_fit(self, x_axis, data, estimator, units=None, add_params=None, **kwargs): +def make_decayexponential_fit(self, x_axis, data, estimator, units=None, + add_params=None, **kwargs): """ Performes a exponential decay with offset fit on the provided data. @param numpy.array x_axis: 1D axis values @@ -224,24 +226,20 @@ def make_decayexponential_fit(self, x_axis, data, estimator, units=None, add_par except: result = exponentialdecay.fit(data, x=x_axis, params=params, **kwargs) self.log.warning('The exponentialdecay with offset fit did not work. ' - 'Message: {}'.format(str(result.message))) + f'Message: {result.message!s}') if units is None: units = ['arb. unit', 'arb. unit'] - result_str_dict = dict() # create result string for gui - - result_str_dict['Amplitude'] = {'value': result.params['amplitude'].value, - 'error': result.params['amplitude'].stderr, - 'unit': units[1]} # amplitude - - result_str_dict['Lifetime'] = {'value': result.params['lifetime'].value, - 'error': result.params['lifetime'].stderr, - 'unit': units[0]} # lifetime - - result_str_dict['Offset'] = {'value': result.params['offset'].value, - 'error': result.params['offset'].stderr, - 'unit': units[1]} # offset + result_str_dict = {'Amplitude': {'value': result.params['amplitude'].value, + 'error': result.params['amplitude'].stderr, + 'unit': units[1]}, + 'Lifetime': {'value': result.params['lifetime'].value, + 'error': result.params['lifetime'].stderr, + 'unit': units[0]}, + 'Offset': {'value': result.params['offset'].value, + 'error': result.params['offset'].stderr, + 'unit': units[1]}} # create result string for gui result.result_str_dict = result_str_dict @@ -270,10 +268,7 @@ def estimate_decayexponential(self, x_axis, data, params): offset = data[-max(1, int(len(x_axis) / 10)):].mean() # substraction of offset, check whether - if data[0] < data[-1]: - data_level = offset - data - else: - data_level = data - offset + data_level = offset - data if data[0] < data[-1] else data - offset # check if the data level contain still negative values and correct # the data level therefore. Otherwise problems in the logarithm appear. @@ -283,7 +278,7 @@ def estimate_decayexponential(self, x_axis, data, params): # remove all the data that can be smaller than or equals to std. # when the data is smaller than std, it is beyond resolution # which is not helpful to our fitting. - for i in range(0, len(x_axis)): + for i in range(len(x_axis)): if data_level[i] <= data_level.std(): break @@ -295,16 +290,21 @@ def estimate_decayexponential(self, x_axis, data, params): data_level_log = np.log(data_level[0:i]) # linear fit, see linearmethods.py - linear_result = self.make_linear_fit(x_axis=x_axis[0:i], data=data_level_log, estimator=self.estimate_linear) - params['lifetime'].set(value=-1 / linear_result.params['slope'].value, min=min_lifetime) + linear_result = self.make_linear_fit(x_axis=x_axis[0:i], data=data_level_log, + estimator=self.estimate_linear) + params['lifetime'].set(value=-1 / linear_result.params['slope'].value, + min=min_lifetime) # amplitude can be positive of negative if data[0] < data[-1]: - params['amplitude'].set(value=-np.exp(linear_result.params['offset'].value), max=-ampl) + params['amplitude'].set(value=-np.exp(linear_result.params['offset'].value), + max=-ampl) else: - params['amplitude'].set(value=np.exp(linear_result.params['offset'].value), min=ampl) + params['amplitude'].set(value=np.exp(linear_result.params['offset'].value), + min=ampl) except: - self.log.warning('Lifetime too small in estimate_exponentialdecay, beyond resolution!') + self.log.warning( + 'Lifetime too small in estimate_exponentialdecay, beyond resolution!') params['lifetime'].set(value=x_axis[i] - x_axis[0], min=min_lifetime) params['amplitude'].set(value=data_level[0]) @@ -318,7 +318,8 @@ def estimate_decayexponential(self, x_axis, data, params): # stretched exponential decay with offset # ############################################# -def make_decayexponentialstretched_fit(self, x_axis, data, estimator, units=None, add_params=None, **kwargs): +def make_decayexponentialstretched_fit(self, x_axis, data, estimator, units=None, + add_params=None, **kwargs): """ Performes a stretched exponential decay with offset fit on the provided data. @param numpy.array x_axis: 1D axis values @@ -345,28 +346,23 @@ def make_decayexponentialstretched_fit(self, x_axis, data, estimator, units=None except: result = stret_exp_decay_offset.fit(data, x=x_axis, params=params, **kwargs) self.log.warning('The double exponentialdecay with offset fit did not work. ' - 'Message: {}'.format(str(result.message))) + f'Message: {result.message!s}') if units is None: units = ['arb. unit', 'arb. unit'] - result_str_dict = dict() # create result string for gui - - result_str_dict['Amplitude'] = {'value': result.params['amplitude'].value, - 'error': result.params['amplitude'].stderr, - 'unit': units[1]} # amplitude - - result_str_dict['Lifetime'] = {'value': result.params['lifetime'].value, - 'error': result.params['lifetime'].stderr, - 'unit': units[0]} # lifetime - - result_str_dict['Offset'] = {'value': result.params['offset'].value, - 'error': result.params['offset'].stderr, - 'unit': units[1]} # offset - - result_str_dict['Beta'] = {'value': result.params['beta'].value, - 'error': result.params['beta'].stderr, - 'unit': ''} # Beta (exponent of exponential exponent) + result_str_dict = {'Amplitude': {'value': result.params['amplitude'].value, + 'error': result.params['amplitude'].stderr, + 'unit': units[1]}, + 'Lifetime': {'value': result.params['lifetime'].value, + 'error': result.params['lifetime'].stderr, + 'unit': units[0]}, + 'Offset': {'value': result.params['offset'].value, + 'error': result.params['offset'].stderr, + 'unit': units[1]}, + 'Beta': {'value': result.params['beta'].value, + 'error': result.params['beta'].stderr, + 'unit': ''}} # create result string for gui result.result_str_dict = result_str_dict @@ -413,7 +409,7 @@ def estimate_decayexponentialstretched(self, x_axis, data, params): # Take all values up to the standard deviation, the remaining values are # more disturbing the estimation then helping: - for stop_index in range(0, len(x_axis)): + for stop_index in range(len(x_axis)): if data_smoothed[stop_index] <= data_smoothed.std(): break @@ -475,39 +471,30 @@ def make_biexponential_fit(self, x_axis, data, estimator, result = model.fit(data, x=x_axis, params=params, **kwargs) except: result = model.fit(data, x=x_axis, params=params, **kwargs) - self.log.warning('The double gaussian dip fit did not work: {0}'.format( - result.message)) + self.log.warning(f'The double gaussian dip fit did not work: {result.message}') # Write the parameters to allow human-readable output to be generated - result_str_dict = dict() - - result_str_dict['1st amplitude'] = {'value': result.params['e0_amplitude'].value, - 'error': result.params['e0_amplitude'].stderr, - 'unit': units[1]} # amplitude - - result_str_dict['1st lifetime'] = {'value': result.params['e0_lifetime'].value, - 'error': result.params['e0_lifetime'].stderr, - 'unit': units[0]} # lifetime - - result_str_dict['1st beta'] = {'value': result.params['e0_beta'].value, - 'error': result.params['e0_beta'].stderr, - 'unit': ''} # Beta (exponent of exponential exponent) - - result_str_dict['2nd amplitude'] = {'value': result.params['e1_amplitude'].value, - 'error': result.params['e1_amplitude'].stderr, - 'unit': units[1]} # amplitude - - result_str_dict['2nd lifetime'] = {'value': result.params['e1_lifetime'].value, - 'error': result.params['e1_lifetime'].stderr, - 'unit': units[0]} # lifetime - - result_str_dict['2nd beta'] = {'value': result.params['e1_beta'].value, - 'error': result.params['e1_beta'].stderr, - 'unit': ''} # Beta (exponent of exponential exponent) - - result_str_dict['offset'] = {'value': result.params['offset'].value, - 'error': result.params['offset'].stderr, - 'unit': units[1]} # offset + result_str_dict = {'1st amplitude': {'value': result.params['e0_amplitude'].value, + 'error': result.params['e0_amplitude'].stderr, + 'unit': units[1]}, + '1st lifetime': {'value': result.params['e0_lifetime'].value, + 'error': result.params['e0_lifetime'].stderr, + 'unit': units[0]}, + '1st beta': {'value': result.params['e0_beta'].value, + 'error': result.params['e0_beta'].stderr, + 'unit': ''}, + '2nd amplitude': {'value': result.params['e1_amplitude'].value, + 'error': result.params['e1_amplitude'].stderr, + 'unit': units[1]}, + '2nd lifetime': {'value': result.params['e1_lifetime'].value, + 'error': result.params['e1_lifetime'].stderr, + 'unit': units[0]}, + '2nd beta': {'value': result.params['e1_beta'].value, + 'error': result.params['e1_beta'].stderr, + 'unit': ''}, + 'offset': {'value': result.params['offset'].value, + 'error': result.params['offset'].stderr, + 'unit': units[1]}} result.result_str_dict = result_str_dict return result @@ -535,10 +522,7 @@ def estimate_biexponential(self, x_axis, data, params): offset = data[-max(1, int(len(x_axis) / 10)):].mean() # substraction of offset, check whether - if data[0] < data[-1]: - data_level = offset - data - else: - data_level = data - offset + data_level = offset - data if data[0] < data[-1] else data - offset # check if the data level contain still negative values and correct # the data level therefore. Otherwise problems in the logarithm appear. @@ -548,7 +532,7 @@ def estimate_biexponential(self, x_axis, data, params): # remove all the data that can be smaller than or equals to std. # when the data is smaller than std, it is beyond resolution # which is not helpful to our fitting. - for i in range(0, len(x_axis)): + for i in range(len(x_axis)): if data_level[i] <= data_level.std(): break @@ -560,19 +544,27 @@ def estimate_biexponential(self, x_axis, data, params): data_level_log = np.log(data_level[0:i]) # linear fit, see linearmethods.py - linear_result = self.make_linear_fit(x_axis=x_axis[0:i], data=data_level_log, estimator=self.estimate_linear) - params['e0_lifetime'].set(value=-1 / linear_result.params['slope'].value, min=min_lifetime) - params['e1_lifetime'].set(value=-1 / linear_result.params['slope'].value, min=min_lifetime) + linear_result = self.make_linear_fit(x_axis=x_axis[0:i], data=data_level_log, + estimator=self.estimate_linear) + params['e0_lifetime'].set(value=-1 / linear_result.params['slope'].value, + min=min_lifetime) + params['e1_lifetime'].set(value=-1 / linear_result.params['slope'].value, + min=min_lifetime) # amplitude can be positive of negative if data[0] < data[-1]: - params['e0_amplitude'].set(value=-np.exp(linear_result.params['offset'].value), max=-ampl) - params['e1_amplitude'].set(value=-np.exp(linear_result.params['offset'].value), max=-ampl) + params['e0_amplitude'].set( + value=-np.exp(linear_result.params['offset'].value), max=-ampl) + params['e1_amplitude'].set( + value=-np.exp(linear_result.params['offset'].value), max=-ampl) else: - params['e0_amplitude'].set(value=np.exp(linear_result.params['offset'].value), min=ampl) - params['e1_amplitude'].set(value=np.exp(linear_result.params['offset'].value), min=ampl) + params['e0_amplitude'].set( + value=np.exp(linear_result.params['offset'].value), min=ampl) + params['e1_amplitude'].set( + value=np.exp(linear_result.params['offset'].value), min=ampl) except: - self.log.warning('Lifetime too small in estimate_exponential, beyond resolution!') + self.log.warning( + 'Lifetime too small in estimate_exponential, beyond resolution!') params['e0_lifetime'].set(value=x_axis[i] - x_axis[0], min=min_lifetime) params['e1_lifetime'].set(value=x_axis[i] - x_axis[0], min=min_lifetime) diff --git a/qudi_hira_analysis/_fitmethods/gaussianlikemethods.py b/qudi_hira_analysis/_fitmethods/gaussianlikemethods.py index ceab3bb..72cd19b 100644 --- a/qudi_hira_analysis/_fitmethods/gaussianlikemethods.py +++ b/qudi_hira_analysis/_fitmethods/gaussianlikemethods.py @@ -1,4 +1,3 @@ -# -*- coding: utf-8 -*- """ This file contains methods for gaussian-like fitting, these methods @@ -21,11 +20,9 @@ top-level directory of this distribution and at """ - from collections import OrderedDict import numpy as np -from lmfit import Parameters from lmfit.models import Model from scipy.ndimage import filters @@ -78,7 +75,7 @@ def physical_gauss(x, center, sigma): amplitude_model, params = self.make_amplitude_model(prefix=prefix) if not isinstance(prefix, str) and prefix is not None: - self.log.error('The passed prefix <{0}> of type {1} is not a string and' + self.log.error('The passed prefix <{}> of type {} is not a string and' 'cannot be used as a prefix and will be ignored for now.' 'Correct that!'.format(prefix, type(prefix))) gaussian_model = Model(physical_gauss, independent_vars=['x']) @@ -90,13 +87,14 @@ def physical_gauss(x, center, sigma): if prefix is None: prefix = '' - full_gaussian_model.set_param_hint('{0!s}fwhm'.format(prefix), - expr="2.3548200450309493*{0}sigma".format(prefix)) + full_gaussian_model.set_param_hint(f'{prefix!s}fwhm', + expr=f"2.3548200450309493*{prefix}sigma") params = full_gaussian_model.make_params() return full_gaussian_model, params + #################################### # 1D Gaussian model with offset # #################################### @@ -120,13 +118,14 @@ def make_gaussian_model(self, prefix=None): if prefix is None: prefix = '' - gaussian_offset_model.set_param_hint('{0}contrast'.format(prefix), - expr='({0}amplitude/offset)*100'.format(prefix)) + gaussian_offset_model.set_param_hint(f'{prefix}contrast', + expr=f'({prefix}amplitude/offset)*100') params = gaussian_offset_model.make_params() return gaussian_offset_model, params + ###################################################### # 1D Gaussian model with linear (inclined) offset # ###################################################### @@ -152,6 +151,7 @@ def make_gaussianlinearoffset_model(self, prefix=None): return gaussian_linear_offset, params + ########################################## # 1D Multiple Gaussian Model with offset # ########################################## @@ -172,25 +172,27 @@ def make_multiplegaussianoffset_model(self, no_of_functions=1): else: prefix = 'g0_' - multi_gaussian_model, params = self.make_gaussianwithoutoffset_model(prefix=prefix) + multi_gaussian_model, params = self.make_gaussianwithoutoffset_model( + prefix=prefix) constant_model, params = self.make_constant_model() multi_gaussian_model = multi_gaussian_model + constant_model - multi_gaussian_model.set_param_hint('{0}contrast'.format(prefix), - expr='({0}amplitude/offset)*100'.format(prefix)) + multi_gaussian_model.set_param_hint(f'{prefix}contrast', + expr=f'({prefix}amplitude/offset)*100') for ii in range(1, no_of_functions): - - prefix = 'g{0:d}_'.format(ii) - multi_gaussian_model += self.make_gaussianwithoutoffset_model(prefix=prefix)[0] - multi_gaussian_model.set_param_hint('{0}contrast'.format(prefix), - expr='({0}amplitude/offset)*100'.format(prefix)) + prefix = f'g{ii:d}_' + multi_gaussian_model += \ + self.make_gaussianwithoutoffset_model(prefix=prefix)[0] + multi_gaussian_model.set_param_hint(f'{prefix}contrast', + expr=f'({prefix}amplitude/offset)*100') params = multi_gaussian_model.make_params() return multi_gaussian_model, params + ########################################## # 1D Double Gaussian Model with offset # ########################################## @@ -205,6 +207,7 @@ def make_gaussiandouble_model(self): return self.make_multiplegaussianoffset_model(no_of_functions=2) + ########################################## # 1D Triple Gaussian Model with offset # ########################################## @@ -219,6 +222,7 @@ def make_gaussiantriple_model(self): return self.make_multiplegaussianoffset_model(no_of_functions=3) + ##################### # 2D gaussian model # ##################### @@ -279,18 +283,19 @@ def twoDgaussian_function(x, amplitude, center_x, center_y, sigma_x, sigma_y, return g.ravel() if not isinstance(prefix, str) and prefix is not None: - self.log.error('The passed prefix <{0}> of type {1} is not a string and' - 'cannot be used as a prefix and will be ignored for now.' - 'Correct that!'.format(prefix, type(prefix))) + self.log.error('The passed prefix <{}> of type {} is not a string and' + 'cannot be used as a prefix and will be ignored for now.' + 'Correct that!'.format(prefix, type(prefix))) gaussian_2d_model = Model(twoDgaussian_function, independent_vars=['x']) else: gaussian_2d_model = Model(twoDgaussian_function, independent_vars=['x'], - prefix=prefix) + prefix=prefix) params = gaussian_2d_model.make_params() return gaussian_2d_model, params + ################################################################################ # # # Fit functions and their estimators # @@ -301,7 +306,8 @@ def twoDgaussian_function(x, amplitude, center_x, center_y, sigma_x, sigma_y, # 1D Gaussian with flat offset # ################################### -def make_gaussian_fit(self, x_axis, data, estimator, units=None, add_params=None, **kwargs): +def make_gaussian_fit(self, x_axis, data, estimator, units=None, add_params=None, + **kwargs): """ Perform a 1D gaussian peak fit on the provided data. @param numpy.array x_axis: 1D axis values @@ -328,38 +334,35 @@ def make_gaussian_fit(self, x_axis, data, estimator, units=None, add_params=None result = mod_final.fit(data, x=x_axis, params=params, **kwargs) except: self.log.warning('The 1D gaussian peak fit did not work. Error ' - 'message: {0}\n'.format(result.message)) + f'message: {result.message}\n') if units is None: - units = ['arb. unit', 'arb. unit'] + units = ['arb. unit', 'arb. unit'] result_str_dict = OrderedDict() # create result string for gui - #result_str_dict['Amplitude'] = {'value': result.params['amplitude'].value, - # 'error': result.params['amplitude'].stderr, + # result_str_dict['Amplitude'] = {'value': result.params['amplitude'].value, # 'unit': units[1]} #amplitude result_str_dict['Position'] = {'value': result.params['center'].value, - 'error': result.params['center'].stderr, - 'unit': units[0]} #position + 'error': result.params['center'].stderr, + 'unit': units[0]} # position - #result_str_dict['Standard deviation'] = {'value': result.params['sigma'].value, - # 'error': result.params['sigma'].stderr, + # result_str_dict['Standard deviation'] = {'value': result.params['sigma'].value, # 'unit': units[0]} #standart deviation result_str_dict['Linewidth'] = {'value': result.params['fwhm'].value, - 'error': result.params['fwhm'].stderr, - 'unit': units[0]} #FWHM + 'error': result.params['fwhm'].stderr, + 'unit': units[0]} # FWHM result_str_dict['Contrast'] = {'value': result.params['contrast'].value, - 'error': result.params['contrast'].stderr, - 'unit': '%'} #Contrast + 'error': result.params['contrast'].stderr, + 'unit': '%'} # Contrast result.result_str_dict = result_str_dict return result - def estimate_gaussian_peak(self, x_axis, data, params): """ Provides a gaussian offset peak estimator. @@ -406,7 +409,7 @@ def estimate_gaussian_peak(self, x_axis, data, params): # calculating the first moment of the gaussian distribution (which is the # mean value), since it is unreliable if the distribution begins or ends at # the edges of the data (but it helps a lot for standard deviation): - mean_val_calc = np.sum(x_axis*data_smoothed) / np.sum(data_smoothed) + mean_val_calc = np.sum(x_axis * data_smoothed) / np.sum(data_smoothed) params['center'].set(value=x_axis[np.argmax(data_smoothed)], min=center_min, max=center_max) @@ -423,18 +426,18 @@ def estimate_gaussian_peak(self, x_axis, data, params): # mean is then higher, which will decrease eventually the initial value for # sigma. But if the peak values is within the distribution the standard # deviation formula performs even better: - params['sigma'].set(value=np.sqrt(abs(mom2 - mean_val_calc**2)), + params['sigma'].set(value=np.sqrt(abs(mom2 - mean_val_calc ** 2)), min=sigma_min, max=sigma_max) - # params['sigma'].set(value=(x_axis.max() - x_axis.min()) / 3.) # Do not set the maximal amplitude value based on the distribution, since # the fit will fail if the peak is at the edges or beyond the range of the # x values. - params['amplitude'].set(value=data_smoothed.max()-data_smoothed.min(), + params['amplitude'].set(value=data_smoothed.max() - data_smoothed.min(), min=ampl_min) return error, params + def estimate_gaussian_dip(self, x_axis, data, params): """ Provides a gaussian offset dip estimator. @@ -457,9 +460,9 @@ def estimate_gaussian_dip(self, x_axis, data, params): data_negative = data * (-1) error, params_ret = self.estimate_gaussian_peak(x_axis, - data_negative, - params_peak - ) + data_negative, + params_peak + ) params['sigma'] = params_ret['sigma'] params['offset'].set(value=-params_ret['offset']) @@ -470,11 +473,13 @@ def estimate_gaussian_dip(self, x_axis, data, params): return error, params + ############################################## # 1D Gaussian with linear inclined offset # ############################################## -def make_gaussianlinearoffset_fit(self, x_axis, data, estimator, units=None, add_params=None, **kwargs): +def make_gaussianlinearoffset_fit(self, x_axis, data, estimator, units=None, + add_params=None, **kwargs): """ Perform a 1D gaussian peak fit with linear offset on the provided data. @param numpy.array x_axis: 1D axis values @param numpy.array data: 1D data, should have the same dimension as x_axis. @@ -499,40 +504,38 @@ def make_gaussianlinearoffset_fit(self, x_axis, data, estimator, units=None, add result = mod_final.fit(data, x=x_axis, params=params, **kwargs) except: self.log.warning('The 1D gaussian peak fit did not work. Error ' - 'message: {0}\n'.format(result.message)) + f'message: {result.message}\n') if units is None: - units = ['arb. unit', 'arb. unit'] + units = ['arb. unit', 'arb. unit'] result_str_dict = OrderedDict() # create result string for gui - #result_str_dict['Amplitude'] = {'value': result.params['amplitude'].value, - # 'error': result.params['amplitude'].stderr, + # result_str_dict['Amplitude'] = {'value': result.params['amplitude'].value, # 'unit': units[1]} #amplitude result_str_dict['Position'] = {'value': result.params['center'].value, - 'error': result.params['center'].stderr, - 'unit': units[0]} #position + 'error': result.params['center'].stderr, + 'unit': units[0]} # position - #result_str_dict['Standard deviation'] = {'value': result.params['sigma'].value, - # 'error': result.params['sigma'].stderr, + # result_str_dict['Standard deviation'] = {'value': result.params['sigma'].value, # 'unit': units[0]} #standart deviation result_str_dict['Linewidth'] = {'value': result.params['fwhm'].value, - 'error': result.params['fwhm'].stderr, - 'unit': units[0]} #FWHM + 'error': result.params['fwhm'].stderr, + 'unit': units[0]} # FWHM result_str_dict['Contrast'] = {'value': result.params['contrast'].value, - 'error': result.params['contrast'].stderr, - 'unit': '%'} #Contrast + 'error': result.params['contrast'].stderr, + 'unit': '%'} # Contrast - #result_str_dict['Slope'] = {'value': result.params['slope'].value, - # 'error': result.params['slope'].stderr, + # result_str_dict['Slope'] = {'value': result.params['slope'].value, # 'unit': ''} #Slope result.result_str_dict = result_str_dict return result + def estimate_gaussianlinearoffset_peak(self, x_axis, data, params): """ Provides a gauss peak estimator with a linear changing offset. @@ -580,6 +583,7 @@ def estimate_gaussianlinearoffset_peak(self, x_axis, data, params): return error, params + ################################# # Two Gaussian with flat offset # ################################# @@ -627,8 +631,7 @@ def make_gaussiandouble_fit(self, x_axis, data, estimator, result = model.fit(data, x=x_axis, params=params, **kwargs) except: result = model.fit(data, x=x_axis, params=params, **kwargs) - self.log.warning('The double gaussian dip fit did not work: {0}'.format( - result.message)) + self.log.warning(f'The double gaussian dip fit did not work: {result.message}') # Write the parameters to allow human-readable output to be generated result_str_dict = OrderedDict() @@ -662,9 +665,10 @@ def make_gaussiandouble_fit(self, x_axis, data, estimator, result.result_str_dict = result_str_dict return result + def estimate_gaussiandouble_peak(self, x_axis, data, params, - threshold_fraction=0.4, minimal_threshold=0.1, - sigma_threshold_fraction=0.2): + threshold_fraction=0.4, minimal_threshold=0.1, + sigma_threshold_fraction=0.2): """ Provide an estimator for a double gaussian peak fit with the parameters coming from the physical properties of an experiment done in gated counter: - positive peak @@ -687,29 +691,29 @@ def estimate_gaussiandouble_peak(self, x_axis, data, params, error = self._check_1D_input(x_axis=x_axis, data=data, params=params) - mod_lor, params_lor = self.make_multiplelorentzian_model(no_of_functions=2) error, params_lor = self.estimate_lorentziandouble_dip(x_axis=x_axis, - data=-data, - params=params_lor, - threshold_fraction=threshold_fraction, - minimal_threshold=minimal_threshold, - sigma_threshold_fraction=sigma_threshold_fraction) + data=-data, + params=params_lor, + threshold_fraction=threshold_fraction, + minimal_threshold=minimal_threshold, + sigma_threshold_fraction=sigma_threshold_fraction) params['g0_amplitude'].value = -params_lor['l0_amplitude'].value params['g0_center'].value = params_lor['l0_center'].value - params['g0_sigma'].value = params_lor['l0_sigma'].value/(np.sqrt(2*np.log(2))) + params['g0_sigma'].value = params_lor['l0_sigma'].value / (np.sqrt(2 * np.log(2))) params['g1_amplitude'].value = -params_lor['l1_amplitude'].value params['g1_center'].value = params_lor['l1_center'].value - params['g1_sigma'].value = params_lor['l1_sigma'].value/(np.sqrt(2*np.log(2))) + params['g1_sigma'].value = params_lor['l1_sigma'].value / (np.sqrt(2 * np.log(2))) params['offset'].value = -params_lor['offset'].value return error, params + def estimate_gaussiandouble_dip(self, x_axis, data, params, - threshold_fraction=0.4, minimal_threshold=0.1, - sigma_threshold_fraction=0.2): + threshold_fraction=0.4, minimal_threshold=0.1, + sigma_threshold_fraction=0.2): """ Provide an estimator for a double gaussian dip fit with the parameters coming from the physical properties of an experiment done in gated counter: - positive peak @@ -735,22 +739,23 @@ def estimate_gaussiandouble_dip(self, x_axis, data, params, mod_lor, params_lor = self.make_multiplelorentzian_model(no_of_functions=2) error, params_lor = self.estimate_lorentziandouble_dip(x_axis=x_axis, - data=data, - params=params_lor, - threshold_fraction=threshold_fraction, - minimal_threshold=minimal_threshold, - sigma_threshold_fraction=sigma_threshold_fraction) + data=data, + params=params_lor, + threshold_fraction=threshold_fraction, + minimal_threshold=minimal_threshold, + sigma_threshold_fraction=sigma_threshold_fraction) params['g0_amplitude'].value = params_lor['l0_amplitude'].value params['g0_center'].value = params_lor['l0_center'].value - params['g0_sigma'].value = params_lor['l0_sigma'].value/(np.sqrt(2*np.log(2))) + params['g0_sigma'].value = params_lor['l0_sigma'].value / (np.sqrt(2 * np.log(2))) params['g1_amplitude'].value = params_lor['l1_amplitude'].value params['g1_center'].value = params_lor['l1_center'].value - params['g1_sigma'].value = params_lor['l1_sigma'].value/(np.sqrt(2*np.log(2))) + params['g1_sigma'].value = params_lor['l1_sigma'].value / (np.sqrt(2 * np.log(2))) params['offset'].value = params_lor['offset'].value return error, params + ################################### # 2D Gaussian with flat offset # ################################### @@ -758,7 +763,8 @@ def estimate_gaussiandouble_dip(self, x_axis, data, params, # TODO: I think this has an offset, and it should be named so to be consistent with # the 1D functions. -def make_twoDgaussian_fit(self, xy_axes, data, estimator, units=None, add_params=None, **kwargs): +def make_twoDgaussian_fit(self, xy_axes, data, estimator, units=None, add_params=None, + **kwargs): """ This method performes a 2D gaussian fit on the provided data. @param numpy.array xy_axes: 2D axes values. xy_axes[0] contains x_axis and @@ -790,11 +796,11 @@ def make_twoDgaussian_fit(self, xy_axes, data, estimator, units=None, add_params result = gaussian_2d_model.fit(data, x=xy_axes, params=params, **kwargs) except: result = gaussian_2d_model.fit(data, x=xy_axes, params=params, **kwargs) - self.log.warning('The 2D gaussian fit did not work: {0}'.format( - result.message)) + self.log.warning(f'The 2D gaussian fit did not work: {result.message}') return result + def estimate_twoDgaussian(self, x_axis, y_axis, data, params): """ Provide a simple two dimensional gaussian function. @@ -815,8 +821,6 @@ def estimate_twoDgaussian(self, x_axis, y_axis, data, params): # FIXME: 1D array x_axis, y_axis, 2D data??? # #needed me 1 hour to think about, but not needed in the end...maybe needed at a later point - # len_x=np.where(x_axis[0]==x_axis)[0][1] - # len_y=len(data)/len_x amplitude = float(data.max() - data.min()) @@ -847,26 +851,27 @@ def estimate_twoDgaussian(self, x_axis, y_axis, data, params): error = -1 # auxiliary variables: - stepsize_x = x_axis[1]-x_axis[0] - stepsize_y = y_axis[1]-y_axis[0] + stepsize_x = x_axis[1] - x_axis[0] + stepsize_y = y_axis[1] - y_axis[0] n_steps_x = len(x_axis) n_steps_y = len(y_axis) # populate the parameter container: params['amplitude'].set(value=amplitude, min=100, max=1e7) - params['sigma_x'].set(value=sigma_x, min=1*stepsize_x, - max=3*(x_axis[-1]-x_axis[0])) - params['sigma_y'].set(value=sigma_y, min=1*stepsize_y, - max=3*(y_axis[-1]-y_axis[0])) - params['center_x'].set(value=center_x, min=(x_axis[0])-n_steps_x*stepsize_x, - max=x_axis[-1]+n_steps_x*stepsize_x) - params['center_y'].set(value=center_y, min=(y_axis[0])-n_steps_y*stepsize_y, - max=y_axis[-1]+n_steps_y*stepsize_y) + params['sigma_x'].set(value=sigma_x, min=1 * stepsize_x, + max=3 * (x_axis[-1] - x_axis[0])) + params['sigma_y'].set(value=sigma_y, min=1 * stepsize_y, + max=3 * (y_axis[-1] - y_axis[0])) + params['center_x'].set(value=center_x, min=(x_axis[0]) - n_steps_x * stepsize_x, + max=x_axis[-1] + n_steps_x * stepsize_x) + params['center_y'].set(value=center_y, min=(y_axis[0]) - n_steps_y * stepsize_y, + max=y_axis[-1] + n_steps_y * stepsize_y) params['theta'].set(value=theta, min=0, max=np.pi) params['offset'].set(value=offset, min=0, max=1e7) return error, params + def estimate_twoDgaussian_MLE(self, x_axis, y_axis, data, params): """ Provide an estimator for 2D gaussian based on maximum likelihood estimation. @@ -919,21 +924,21 @@ def estimate_twoDgaussian_MLE(self, x_axis, y_axis, data, params): error = -1 # auxiliary variables: - stepsize_x = x_axis[1]-x_axis[0] - stepsize_y = y_axis[1]-y_axis[0] + stepsize_x = x_axis[1] - x_axis[0] + stepsize_y = y_axis[1] - y_axis[0] n_steps_x = len(x_axis) n_steps_y = len(y_axis) # populate the parameter container: params['amplitude'].set(value=amplitude, min=100, max=1e7) - params['sigma_x'].set(value=sigma_x, min=1*stepsize_x, - max=3*(x_axis[-1]-x_axis[0])) - params['sigma_y'].set(value=sigma_y, min=1*stepsize_y, - max=3*(y_axis[-1]-y_axis[0])) - params['center_x'].set(value=center_x, min=(x_axis[0])-n_steps_x*stepsize_x, - max=x_axis[-1]+n_steps_x*stepsize_x) - params['center_y'].set(value=center_y, min=(y_axis[0])-n_steps_y*stepsize_y, - max=y_axis[-1]+n_steps_y*stepsize_y) + params['sigma_x'].set(value=sigma_x, min=1 * stepsize_x, + max=3 * (x_axis[-1] - x_axis[0])) + params['sigma_y'].set(value=sigma_y, min=1 * stepsize_y, + max=3 * (y_axis[-1] - y_axis[0])) + params['center_x'].set(value=center_x, min=(x_axis[0]) - n_steps_x * stepsize_x, + max=x_axis[-1] + n_steps_x * stepsize_x) + params['center_y'].set(value=center_y, min=(y_axis[0]) - n_steps_y * stepsize_y, + max=y_axis[-1] + n_steps_y * stepsize_y) params['theta'].set(value=theta, min=0, max=np.pi) params['offset'].set(value=offset, min=0, max=1e7) diff --git a/qudi_hira_analysis/_fitmethods/generalmethods.py b/qudi_hira_analysis/_fitmethods/generalmethods.py index db82d52..d1422d5 100644 --- a/qudi_hira_analysis/_fitmethods/generalmethods.py +++ b/qudi_hira_analysis/_fitmethods/generalmethods.py @@ -1,4 +1,3 @@ -# -*- coding: utf-8 -*- """ This file contains general methods which are imported by class FitLogic. These general methods are available for all different fitting methods. @@ -20,7 +19,6 @@ top-level directory of this distribution and at """ - from collections import OrderedDict import lmfit @@ -87,11 +85,13 @@ def _substitute_params(self, initial_params, update_params=None): if update_params[para].value is not None: # Adapt the limits to the value: - if (initial_params[para].min is not None) and (initial_params[para].min > update_params[para].value): - initial_params[para].min = update_params[para].value + if (initial_params[para].min is not None) and ( + initial_params[para].min > update_params[para].value): + initial_params[para].min = update_params[para].value - if (initial_params[para].max is not None) and (initial_params[para].max < update_params[para].value): - initial_params[para].max = update_params[para].value + if (initial_params[para].max is not None) and ( + initial_params[para].max < update_params[para].value): + initial_params[para].max = update_params[para].value initial_params[para].value = update_params[para].value @@ -116,22 +116,25 @@ def _substitute_params(self, initial_params, update_params=None): if 'value' in update_params[para]: # Adapt the limits to the value: - if (initial_params[para].min is not None) and (initial_params[para].min > update_params[para]['value']): + if (initial_params[para].min is not None) and ( + initial_params[para].min > update_params[para]['value']): initial_params[para].min = update_params[para]['value'] - if (initial_params[para].max is not None) and (initial_params[para].max < update_params[para]['value']): + if (initial_params[para].max is not None) and ( + initial_params[para].max < update_params[para]['value']): initial_params[para].max = update_params[para]['value'] initial_params[para].value = update_params[para]['value'] else: - self.log.error('The type of the passed update_params object <{0}> is ' - 'neither of type lmfit.parameter.Parameters, ' - 'OrderedDict or dict! Correct that, the initial_params' - 'will be returned.'.format(type(update_params))) + self.log.error('The type of the passed update_params object <{}> is ' + 'neither of type lmfit.parameter.Parameters, ' + 'OrderedDict or dict! Correct that, the initial_params' + 'will be returned.'.format(type(update_params))) return initial_params + def create_fit_string(self, result, model, units=None, decimal_digits_value_given=None, decimal_digits_err_given=None): """ This method can produces a well readable string from the results of a fitted model. @@ -148,7 +151,7 @@ def create_fit_string(self, result, model, units=None, decimal_digits_value_give @return str fit_result: readable string """ if units is None: - units = dict() + units = {} # TODO: Add multiplicator # TODO: Add decimal dict # TODO: Add sensible output such that e only multiple of 3 and err and value have same exponent @@ -156,10 +159,10 @@ def create_fit_string(self, result, model, units=None, decimal_digits_value_give fit_result = '' for variable in model.param_names: # check order of number - exponent_error = int("{:e}".format(result.params[variable].stderr)[-3:]) - exponent_value = int("{:e}".format(result.params[variable].value)[-3:]) + exponent_error = int(f"{result.params[variable].stderr:e}"[-3:]) + exponent_value = int(f"{result.params[variable].value:e}"[-3:]) if decimal_digits_value_given is None: - decimal_digits_value = int(exponent_value-exponent_error+1) + decimal_digits_value = int(exponent_value - exponent_error + 1) if decimal_digits_value <= 0: decimal_digits_value = 1 if decimal_digits_err_given is None: @@ -167,87 +170,93 @@ def create_fit_string(self, result, model, units=None, decimal_digits_value_give if decimal_digits_err <= 0: decimal_digits_err = 1 try: - fit_result += ("{0} [{1}] : {2} ± {3}\n".format(str(variable), - units[variable], - "{0:.{1}e}".format( - float(result.params[variable].value), - decimal_digits_value), - "{0:.{1}e}".format( - float(result.params[variable].stderr), - decimal_digits_err))) + fit_result += ("{} [{}] : {} ± {}\n".format(str(variable), + units[variable], + "{0:.{1}e}".format( + float(result.params[ + variable].value), + decimal_digits_value), + "{0:.{1}e}".format( + float(result.params[ + variable].stderr), + decimal_digits_err))) except: # self.log.warning('No unit given for parameter {}, setting unit ' # 'to empty string'.format(variable)) - fit_result += ("{0} [{1}] : {2} ± {3}\n".format(str(variable), - "arb. u.", - "{0:.{1}e}".format( - float(result.params[variable].value), - decimal_digits_value), - "{0:.{1}e}".format( - float(result.params[variable].stderr), - decimal_digits_err))) + fit_result += ("{} [{}] : {} ± {}\n".format(str(variable), + "arb. u.", + "{0:.{1}e}".format( + float(result.params[ + variable].value), + decimal_digits_value), + "{0:.{1}e}".format( + float(result.params[ + variable].stderr), + decimal_digits_err))) return fit_result -def _search_end_of_dip(self, direction, data, peak_arg, start_arg, end_arg, sigma_threshold, minimal_threshold, make_prints): +def _search_end_of_dip(self, direction, data, peak_arg, start_arg, end_arg, + sigma_threshold, + minimal_threshold, make_prints): """ data has to be offset leveled such that offset is substracted """ # Todo: Create doc string - absolute_min = data[peak_arg] + absolute_min = data[peak_arg] if direction == 'left': mult = -1 - sigma_arg=start_arg + sigma_arg = start_arg elif direction == 'right': mult = +1 - sigma_arg=end_arg + sigma_arg = end_arg else: print('No valid direction in search end of peak') - ii=0 + ii = 0 - #if the minimum is at the end set this as boarder - if (peak_arg != start_arg and direction=='left' or - peak_arg != end_arg and direction=='right'): + # if the minimum is at the end set this as boarder + if (peak_arg != start_arg and direction == 'left' or + peak_arg != end_arg and direction == 'right'): while True: # if no minimum can be found decrease threshold - if ((peak_arg-iiend_arg and direction=='right')): - sigma_threshold*=0.9 - ii=0 + if ((peak_arg - ii < start_arg and direction == 'left') or + (peak_arg + ii > end_arg and direction == 'right')): + sigma_threshold *= 0.9 + ii = 0 if make_prints: - print('h1 sigma_threshold',sigma_threshold) + print('h1 sigma_threshold', sigma_threshold) - #if the dip is always over threshold the end is as + # if the dip is always over threshold the end is as # set before - if abs(sigma_threshold/absolute_min) abs(threshold) and \ - abs(left_min) > abs(right_min): + abs(left_min) > abs(right_min): if make_prints: print('h13') # there is a minimum on the left side which is higher # than the minimum on the right side - dip1_arg = left_argmin+left_index + dip1_arg = left_argmin + left_index break - elif abs(right_min)>abs(threshold): + elif abs(right_min) > abs(threshold): # there is a minimum on the right side which is higher # than on left side - dip1_arg=right_argmin+mid_index_right + dip1_arg = right_argmin + mid_index_right if make_prints: print('h14') break else: # no minimum at all over threshold so lowering threshold # and resetting search area - threshold*=0.9 - left_index=int(0) - right_index=len(x_axis)-1 - mid_index_left=sigma0_argleft - mid_index_right=sigma0_argright + threshold *= 0.9 + left_index = 0 + right_index = len(x_axis) - 1 + mid_index_left = sigma0_argleft + mid_index_right = sigma0_argright if make_prints: print('h15') - #if no second dip can be found set both to same value - if abs(threshold/absolute_min) distance_right: - dip1_arg = dip0_arg - abs(distance_left-distance_right) + dip1_arg = dip0_arg - abs(distance_left - distance_right) elif distance_left < distance_right: - dip1_arg = dip0_arg + abs(distance_left-distance_right) + dip1_arg = dip0_arg + abs(distance_left - distance_right) else: dip1_arg = dip0_arg - #print(distance_left,distance_right,dip1_arg) else: # if the peaks are not overlapping search for left and right # boarder of the dip # ====== search for the right end of the dip ====== sigma_threshold, sigma1_argleft = self._search_end_of_dip( - direction='left', - data=data, - peak_arg = dip1_arg, - start_arg = 0, - end_arg = len(data)-1, - sigma_threshold = sigma_threshold_fraction*absolute_min, - minimal_threshold = minimal_threshold, - make_prints= make_prints) + direction='left', + data=data, + peak_arg=dip1_arg, + start_arg=0, + end_arg=len(data) - 1, + sigma_threshold=sigma_threshold_fraction * absolute_min, + minimal_threshold=minimal_threshold, + make_prints=make_prints) # ====== search for the right end of the dip ====== sigma_threshold, sigma1_argright = self._search_end_of_dip( - direction='right', - data=data, - peak_arg = dip1_arg, - start_arg = 0, - end_arg = len(data)-1, - sigma_threshold = sigma_threshold_fraction*absolute_min, - minimal_threshold = minimal_threshold, - make_prints= make_prints) + direction='right', + data=data, + peak_arg=dip1_arg, + start_arg=0, + end_arg=len(data) - 1, + sigma_threshold=sigma_threshold_fraction * absolute_min, + minimal_threshold=minimal_threshold, + make_prints=make_prints) return error, sigma0_argleft, dip0_arg, sigma0_argright, sigma1_argleft, dip1_arg, sigma1_argright @@ -471,7 +478,7 @@ def find_offset_parameter(self, x_values=None, data=None): elif len(x_values) >= 100.: len_x = 10 else: - len_x = int(len(x_values)/10.)+1 + len_x = int(len(x_values) / 10.) + 1 lorentz = mod.eval(x=np.linspace(0, len_x, len_x), amplitude=1, offset=0., sigma=len_x / 4., center=len_x / 2.) @@ -480,10 +487,11 @@ def find_offset_parameter(self, x_values=None, data=None): # finding most frequent value which is supposed to be the offset hist = np.histogram(data_smooth, bins=10) - offset = (hist[1][hist[0].argmax()]+hist[1][hist[0].argmax()+1])/2. + offset = (hist[1][hist[0].argmax()] + hist[1][hist[0].argmax() + 1]) / 2. return data_smooth, offset + ############################################################################ # # # Additional routines with gaussian-like filter # @@ -501,7 +509,7 @@ def gaussian_smoothing(self, data=None, filter_len=None, filter_sigma=None): @return array: smoothed data """ - #Todo: Check for wrong data type + # Todo: Check for wrong data type if filter_len is None: if len(data) < 20.: filter_len = 5 @@ -516,7 +524,6 @@ def gaussian_smoothing(self, data=None, filter_len=None, filter_sigma=None): return convolve1d(data, gaus / gaus.sum(), mode='mirror') - def _check_1D_input(self, x_axis, data, params): """ Helper function to check the input of the fit for general consistency. @@ -542,4 +549,3 @@ def _check_1D_input(self, x_axis, data, params): error = -1 return error - diff --git a/qudi_hira_analysis/_fitmethods/hyperbolicsaturationmethods.py b/qudi_hira_analysis/_fitmethods/hyperbolicsaturationmethods.py index a162e19..3904867 100644 --- a/qudi_hira_analysis/_fitmethods/hyperbolicsaturationmethods.py +++ b/qudi_hira_analysis/_fitmethods/hyperbolicsaturationmethods.py @@ -1,4 +1,3 @@ -# -*- coding: utf-8 -*- """ This file contains methods for hyperbolic saturation fitting, these methods are imported by class FitLogic. @@ -20,11 +19,9 @@ top-level directory of this distribution and at """ - import numpy as np from lmfit.models import Model - ################################################################################ # # # Hyperbolic saturation models # @@ -63,9 +60,9 @@ def hyperbolicsaturation_function(x, I_sat, P_sat): return I_sat * (x / (x + P_sat)) if not isinstance(prefix, str) and prefix is not None: - self.log.error('The passed prefix <{0}> of type {1} is not a string and' - 'cannot be used as a prefix and will be ignored for now.' - 'Correct that!'.format(prefix, type(prefix))) + self.log.error('The passed prefix <{}> of type {} is not a string and' + 'cannot be used as a prefix and will be ignored for now.' + 'Correct that!'.format(prefix, type(prefix))) mod_sat = Model(hyperbolicsaturation_function, independent_vars=['x']) else: @@ -80,7 +77,8 @@ def hyperbolicsaturation_function(x, I_sat, P_sat): return complete_model, params -def make_hyperbolicsaturation_fit(self, x_axis, data, estimator, units=None, add_params=None, **kwargs): +def make_hyperbolicsaturation_fit(self, x_axis, data, estimator, units=None, + add_params=None, **kwargs): """ Perform a fit on the provided data with a fluorescence depending function. @param numpy.array x_axis: 1D axis values @@ -128,23 +126,22 @@ def estimate_hyperbolicsaturation(self, x_axis, data, params): error = self._check_1D_input(x_axis=x_axis, data=data, params=params) - x_axis_half = x_axis[len(x_axis)//2:] - data_half = data[len(x_axis)//2:] + x_axis_half = x_axis[len(x_axis) // 2:] + data_half = data[len(x_axis) // 2:] results_lin = self.make_linear_fit(x_axis=x_axis_half, data=data_half, - estimator=self.estimate_linear) + estimator=self.estimate_linear) est_slope = results_lin.params['slope'].value est_offset = data.min() - data_red = data - est_slope*x_axis - est_offset - est_I_sat = np.mean(data_red[len(data_red)//2:]) - est_P_sat = est_I_sat/2 + data_red = data - est_slope * x_axis - est_offset + est_I_sat = np.mean(data_red[len(data_red) // 2:]) + est_P_sat = est_I_sat / 2 params['I_sat'].value = est_I_sat params['slope'].value = est_slope params['offset'].value = est_offset params['P_sat'].value = est_P_sat - return error, params diff --git a/qudi_hira_analysis/_fitmethods/linearmethods.py b/qudi_hira_analysis/_fitmethods/linearmethods.py index 0a60902..aaf64d9 100644 --- a/qudi_hira_analysis/_fitmethods/linearmethods.py +++ b/qudi_hira_analysis/_fitmethods/linearmethods.py @@ -1,4 +1,3 @@ -# -*- coding: utf-8 -*- """ This file contains methods for linear fitting, these methods are imported by class FitLogic. The functions can be used for amy @@ -55,6 +54,7 @@ def make_constant_model(self, prefix=None): For further information have a look in: http://cars9.uchicago.edu/software/python/lmfit/builtin_models.html#models.GaussianModel """ + def constant_function(x, offset): """ Function of a constant value. @@ -67,9 +67,10 @@ def constant_function(x, offset): return offset if not isinstance(prefix, str) and prefix is not None: - self.log.error('The passed prefix <{0}> of type {1} is not a string and cannot be used as ' - 'a prefix and will be ignored for now. Correct that!'.format(prefix, - type(prefix))) + self.log.error( + 'The passed prefix <{}> of type {} is not a string and cannot be used as ' + 'a prefix and will be ignored for now. Correct that!'.format(prefix, + type(prefix))) model = Model(constant_function, independent_vars=['x']) else: model = Model(constant_function, independent_vars=['x'], prefix=prefix) @@ -103,9 +104,10 @@ def amplitude_function(x, amplitude): return amplitude if not isinstance(prefix, str) and prefix is not None: - self.log.error('The passed prefix <{0}> of type {1} is not a string and cannot be used as ' - 'a prefix and will be ignored for now. Correct that!'.format(prefix, - type(prefix))) + self.log.error( + 'The passed prefix <{}> of type {} is not a string and cannot be used as ' + 'a prefix and will be ignored for now. Correct that!'.format(prefix, + type(prefix))) model = Model(amplitude_function, independent_vars=['x']) else: model = Model(amplitude_function, independent_vars=['x'], prefix=prefix) @@ -139,9 +141,10 @@ def slope_function(x, slope): return slope if not isinstance(prefix, str) and prefix is not None: - self.log.error('The passed prefix <{0}> of type {1} is not a string and cannot be used as ' - 'a prefix and will be ignored for now. Correct that!'.format(prefix, - type(prefix))) + self.log.error( + 'The passed prefix <{}> of type {} is not a string and cannot be used as ' + 'a prefix and will be ignored for now. Correct that!'.format(prefix, + type(prefix))) model = Model(slope_function, independent_vars=['x']) else: model = Model(slope_function, independent_vars=['x'], prefix=prefix) @@ -174,9 +177,10 @@ def linear_function(x): return x if not isinstance(prefix, str) and prefix is not None: - self.log.error('The passed prefix <{0}> of type {1} is not a string and cannot be used as ' - 'a prefix and will be ignored for now. Correct that!'.format(prefix, - type(prefix))) + self.log.error( + 'The passed prefix <{}> of type {} is not a string and cannot be used as ' + 'a prefix and will be ignored for now. Correct that!'.format(prefix, + type(prefix))) linear_mod = Model(linear_function, independent_vars=['x']) else: linear_mod = Model(linear_function, independent_vars=['x'], prefix=prefix) @@ -190,7 +194,8 @@ def linear_function(x): return model, params -def make_linear_fit(self, x_axis, data, estimator, units=None, add_params=None, **kwargs): +def make_linear_fit(self, x_axis, data, estimator, units=None, add_params=None, + **kwargs): """ Performe a linear fit on the provided data. @param numpy.array x_axis: 1D axis values @@ -218,11 +223,11 @@ def make_linear_fit(self, x_axis, data, estimator, units=None, add_params=None, if units is None: units = ['arb. unit', 'arb. unit'] - result_str_dict = dict() + result_str_dict = {} result_str_dict['Slope'] = {'value': result.params['slope'].value, 'error': result.params['slope'].stderr, - 'unit': '{0}/{1}'.format(units[1], units[0])} + 'unit': f'{units[1]}/{units[0]}'} result_str_dict['Offset'] = {'value': result.params['offset'].value, 'error': result.params['offset'].stderr, 'unit': units[1]} @@ -256,10 +261,10 @@ def estimate_linear(self, x_axis, data, params): x_mean = x_axis.mean() data_mean = data.mean() - for i in range(0, len(x_axis)): - a_1 += (x_axis[i]-x_mean)*(data[i]-data_mean) - a_2 += np.power(x_axis[i]-x_mean, 2) - slope = a_1/a_2 + for i in range(len(x_axis)): + a_1 += (x_axis[i] - x_mean) * (data[i] - data_mean) + a_2 += np.power(x_axis[i] - x_mean, 2) + slope = a_1 / a_2 intercept = data_mean - slope * x_mean params['offset'].value = intercept params['slope'].value = slope @@ -269,4 +274,3 @@ def estimate_linear(self, x_axis, data, params): params['offset'].value = 0 return error, params - diff --git a/qudi_hira_analysis/_fitmethods/lorentzianlikemethods.py b/qudi_hira_analysis/_fitmethods/lorentzianlikemethods.py index 82a9aa2..f1eb5e3 100644 --- a/qudi_hira_analysis/_fitmethods/lorentzianlikemethods.py +++ b/qudi_hira_analysis/_fitmethods/lorentzianlikemethods.py @@ -1,4 +1,3 @@ -# -*- coding: utf-8 -*- """ This file contains methods for lorentzian-like fitting, these methods are imported by class FitLogic. @@ -22,11 +21,9 @@ top-level directory of this distribution and at """ - from collections import OrderedDict import numpy as np -from lmfit import Parameters from lmfit.models import Model from scipy.interpolate import InterpolatedUnivariateSpline from scipy.ndimage import filters @@ -100,6 +97,7 @@ """ + #################################### # Lorentzian model # #################################### @@ -145,9 +143,9 @@ def physical_lorentzian(x, center, sigma): if not isinstance(prefix, str) and prefix is not None: self.log.error( - 'The passed prefix <{0}> of type {1} is not a string and' + f'The passed prefix <{prefix}> of type {type(prefix)} is not a string and' 'cannot be used as a prefix and will be ignored for now.' - 'Correct that!'.format(prefix, type(prefix))) + 'Correct that!') lorentz_model = Model(physical_lorentzian, independent_vars=['x']) else: lorentz_model = Model( @@ -163,12 +161,11 @@ def physical_lorentzian(x, center, sigma): if prefix is None: prefix = '' full_lorentz_model.set_param_hint( - '{0!s}fwhm'.format(prefix), - expr="2*{0!s}sigma".format(prefix)) + f'{prefix!s}fwhm', + expr=f"2*{prefix!s}sigma") # full_lorentz_model.set_param_hint('{0}contrast'.format(prefix), # expr='(-100.0)') - # expr='({0!s}amplitude/offset)*100'.format(prefix)) - # params.add('{0}contrast'.format(prefix), expr='({0!s}amplitude/offset)*100'.format(prefix)) + # expr='({0!s}amplitude/offset)*100'.format(prefix)) return full_lorentz_model, params @@ -197,8 +194,8 @@ def make_lorentzian_model(self, prefix=None): if prefix is None: prefix = '' - lorentz_offset_model.set_param_hint('{0}contrast'.format(prefix), - expr='({0}amplitude/offset)*100'.format(prefix)) + lorentz_offset_model.set_param_hint(f'{prefix}contrast', + expr=f'({prefix}amplitude/offset)*100') params = lorentz_offset_model.make_params() @@ -223,27 +220,29 @@ def make_multiplelorentzian_model(self, no_of_functions=1): multi_lorentz_model, params = self.make_lorentzian_model() else: prefix = 'l0_' - multi_lorentz_model, params = self.make_lorentzianwithoutoffset_model(prefix=prefix) + multi_lorentz_model, params = self.make_lorentzianwithoutoffset_model( + prefix=prefix) constant_model, params = self.make_constant_model() multi_lorentz_model = multi_lorentz_model + constant_model multi_lorentz_model.set_param_hint( - '{0}contrast'.format(prefix), - expr='({0}amplitude/offset)*100'.format(prefix)) - + f'{prefix}contrast', + expr=f'({prefix}amplitude/offset)*100') for ii in range(1, no_of_functions): - prefix = 'l{0:d}_'.format(ii) - multi_lorentz_model += self.make_lorentzianwithoutoffset_model(prefix=prefix)[0] + prefix = f'l{ii:d}_' + multi_lorentz_model += \ + self.make_lorentzianwithoutoffset_model(prefix=prefix)[0] multi_lorentz_model.set_param_hint( - '{0}contrast'.format(prefix), - expr='({0}amplitude/offset)*100'.format(prefix)) + f'{prefix}contrast', + expr=f'({prefix}amplitude/offset)*100') params = multi_lorentz_model.make_params() return multi_lorentz_model, params + ################################################# # Double Lorentzian model with offset # ################################################# @@ -257,6 +256,7 @@ def make_lorentziandouble_model(self): return self.make_multiplelorentzian_model(no_of_functions=2) + ################################################# # Triple Lorentzian model with offset # ################################################# @@ -270,6 +270,7 @@ def make_lorentziantriple_model(self): return self.make_multiplelorentzian_model(no_of_functions=3) + ################################################################################ # # # Fit functions and their estimators # @@ -309,7 +310,7 @@ def make_lorentzian_fit(self, x_axis, data, estimator, units=None, except: result = model.fit(data, x=x_axis, params=params, **kwargs) self.log.warning('The 1D lorentzian fit did not work. Error ' - 'message: {0}\n'.format(result.message)) + f'message: {result.message}\n') # Write the parameters to allow human-readable output to be generated result_str_dict = OrderedDict() @@ -334,6 +335,7 @@ def make_lorentzian_fit(self, x_axis, data, estimator, units=None, result.result_str_dict = result_str_dict return result + def estimate_lorentzian_dip(self, x_axis, data, params): """ Provides an estimator to obtain initial values for lorentzian function. @@ -359,15 +361,14 @@ def estimate_lorentzian_dip(self, x_axis, data, params): data_smooth, offset = self.find_offset_parameter(x_axis, data) - # data_level = data-offset data_level = data_smooth - offset # calculate from the leveled data the amplitude: amplitude = data_level.min() - smoothing_spline = 1 # must be 1<= smoothing_spline <= 5 + smoothing_spline = 1 # must be 1<= smoothing_spline <= 5 fit_function = InterpolatedUnivariateSpline(x_axis, data_level, - k=smoothing_spline) + k=smoothing_spline) numerical_integral = fit_function.integral(x_axis[0], x_axis[-1]) x_zero = x_axis[np.argmin(data_smooth)] @@ -390,7 +391,8 @@ def estimate_lorentzian_dip(self, x_axis, data, params): return error, params -def estimate_lorentzian_peak (self, x_axis, data, params): + +def estimate_lorentzian_peak(self, x_axis, data, params): """ Provides a lorentzian offset peak estimator. @param numpy.array x_axis: 1D axis values @@ -433,7 +435,8 @@ def estimate_lorentzian_peak (self, x_axis, data, params): # Double Lorentzian with offset fitting # ################################################################################ -def make_lorentziandouble_fit(self, x_axis, data, estimator, units=None, add_params=None, **kwargs): +def make_lorentziandouble_fit(self, x_axis, data, estimator, units=None, + add_params=None, **kwargs): """ Perform a 1D double lorentzian dip fit with offset on the provided data. @param numpy.array x_axis: 1D axis values @@ -463,7 +466,7 @@ def make_lorentziandouble_fit(self, x_axis, data, estimator, units=None, add_par except: result = model.fit(data, x=x_axis, params=params, **kwargs) self.log.error('The double lorentzian fit did not ' - 'work: {0}'.format(result.message)) + f'work: {result.message}') # Write the parameters to allow human-readable output to be generated result_str_dict = OrderedDict() @@ -512,6 +515,7 @@ def make_lorentziandouble_fit(self, x_axis, data, estimator, units=None, add_par result.result_str_dict = result_str_dict return result + def estimate_lorentziandouble_dip(self, x_axis, data, params, threshold_fraction=0.3, minimal_threshold=0.01, @@ -562,11 +566,11 @@ def estimate_lorentziandouble_dip(self, x_axis, data, params, # (x_axis[sigma0_argright] - x_axis[sigma0_argleft]) / # len(data_level[sigma0_argleft:sigma0_argright])) - smoothing_spline = 1 # must be 1<= smoothing_spline <= 5 + smoothing_spline = 1 # must be 1<= smoothing_spline <= 5 fit_function = InterpolatedUnivariateSpline(x_axis, data_level, - k=smoothing_spline) + k=smoothing_spline) numerical_integral_0 = fit_function.integral(x_axis[sigma0_argleft], - x_axis[sigma0_argright]) + x_axis[sigma0_argright]) lorentz0_sigma = abs(numerical_integral_0 / (np.pi * lorentz0_amplitude)) @@ -575,8 +579,6 @@ def estimate_lorentziandouble_dip(self, x_axis, data, params, lorentz1_sigma = abs(numerical_integral_1 / (np.pi * lorentz1_amplitude)) # esstimate amplitude - # lorentz0_amplitude = -1*abs(lorentz0_amplitude*np.pi*lorentz0_sigma) - # lorentz1_amplitude = -1*abs(lorentz1_amplitude*np.pi*lorentz1_sigma) stepsize = x_axis[1] - x_axis[0] full_width = x_axis[-1] - x_axis[0] @@ -613,6 +615,7 @@ def estimate_lorentziandouble_dip(self, x_axis, data, params, return error, params + def estimate_lorentziandouble_peak(self, x_axis, data, params, threshold_fraction=0.3, minimal_threshold=0.01, @@ -684,16 +687,16 @@ def estimate_lorentziandouble_N15(self, x_axis, data, params): # check if parameters make sense error = self._check_1D_input(x_axis=x_axis, data=data, params=params) - hf_splitting = 3.03 * 1e6 # Hz + hf_splitting = 3.03 * 1e6 # Hz # this is an estimator, for a physical application, therefore the x_axis # should fulfill certain constraints: length_x_scan = x_axis[-1] - x_axis[0] - if length_x_scan < hf_splitting/2 or hf_splitting > 1e9: + if length_x_scan < hf_splitting / 2 or hf_splitting > 1e9: self.log.error('The N15 estimator expects an x_axis with a length in the ' - 'range [{0},{1}]Hz, but the passed x_axis has a length of ' - '{2}, which is not sensible for the N15 estimator. Correct ' + 'range [{},{}]Hz, but the passed x_axis has a length of ' + '{}, which is not sensible for the N15 estimator. Correct ' 'that!'.format(hf_splitting / 2, 1e9, length_x_scan)) return -1, params @@ -705,7 +708,8 @@ def estimate_lorentziandouble_N15(self, x_axis, data, params): # filter should have a width of 4 MHz x_filter = np.linspace(0, 4 * points_within_1MHz, 4 * points_within_1MHz) lorentz = np.piecewise(x_filter, [(x_filter >= 0) * (x_filter < len(x_filter) / 4), - (x_filter >= len(x_filter) / 4) * (x_filter < len(x_filter) * 3 / 4), + (x_filter >= len(x_filter) / 4) * ( + x_filter < len(x_filter) * 3 / 4), (x_filter >= len(x_filter) * 3 / 4)], [1, 0, 1]) @@ -719,7 +723,6 @@ def estimate_lorentziandouble_N15(self, x_axis, data, params): else: x_axis_min = x_axis[data_smooth_lorentz.argmin()] - # data_level = data_smooth_lorentz - data_smooth_lorentz.max() data_level = data_smooth_lorentz - offset minimum_level = data_level.min() @@ -743,7 +746,7 @@ def estimate_lorentziandouble_N15(self, x_axis, data, params): params['l1_amplitude'].set(value=params['l0_amplitude'].value, max=-1e-6) params['l1_center'].set(value=params['l0_center'].value + hf_splitting, - expr='l0_center+{0}'.format(hf_splitting)) + expr=f'l0_center+{hf_splitting}') params['l1_sigma'].set(value=params['l0_sigma'].value, min=minimal_sigma, max=maximal_sigma, expr='l0_sigma') @@ -757,12 +760,12 @@ def estimate_lorentziandouble_N15(self, x_axis, data, params): # Triple Lorentzian fitting # # # ############################################################################ -#Todo: check where code breaks +# Todo: check where code breaks # Old Method Names: # make_N14_fit def make_lorentziantriple_fit(self, x_axis, data, estimator, units=None, - add_params=None, **kwargs): + add_params=None, **kwargs): """ Perform a triple lorentzian fit @param numpy.array x_axis: 1D axis values @@ -790,7 +793,7 @@ def make_lorentziantriple_fit(self, x_axis, data, estimator, units=None, except: result = model.fit(data, x=x_axis, params=params, **kwargs) self.log.error('The triple lorentzian fit did not ' - 'work: {0}'.format(result.message)) + f'work: {result.message}') # Write the parameters to allow human-readable output to be generated result_str_dict = OrderedDict() @@ -839,6 +842,7 @@ def make_lorentziantriple_fit(self, x_axis, data, estimator, units=None, result.result_str_dict = result_str_dict return result + def estimate_lorentziantriple_N14(self, x_axis, data, params): """ Estimation of a the hyperfine interaction of a N14 nuclear spin. @@ -868,17 +872,17 @@ def estimate_lorentziantriple_N14(self, x_axis, data, params): # check if parameters make sense error = self._check_1D_input(x_axis=x_axis, data=data, params=params) - hf_splitting = 2.15e6 # hyperfine splitting for a N14 spin + hf_splitting = 2.15e6 # hyperfine splitting for a N14 spin # this is an estimator, for a physical application, therefore the x_axis # should fulfill certain constraints: length_x_scan = x_axis[-1] - x_axis[0] - if length_x_scan < hf_splitting/2 or hf_splitting > 1e9: + if length_x_scan < hf_splitting / 2 or hf_splitting > 1e9: self.log.error('The N14 estimator expects an x_axis with a length in the ' - 'range [{0},{1}]Hz, but the passed x_axis has a length of ' - '{2}, which is not sensible for the N14 estimator. Correct ' - 'that!'.format(hf_splitting/2, 1e9, length_x_scan)) + 'range [{},{}]Hz, but the passed x_axis has a length of ' + '{}, which is not sensible for the N14 estimator. Correct ' + 'that!'.format(hf_splitting / 2, 1e9, length_x_scan)) return -1, params # find the offset parameter, which should be in the fit the zero level: @@ -890,26 +894,30 @@ def estimate_lorentziantriple_N14(self, x_axis, data, params): # filter. Take that to obtain from that the accurate peak position: # filter of one dip should always have a length of approx linewidth 1MHz - points_within_1MHz = len(x_axis)/(x_axis.max()-x_axis.min()) * 1e6 + points_within_1MHz = len(x_axis) / (x_axis.max() - x_axis.min()) * 1e6 # filter should have a width of 5MHz - x_filter = np.linspace(0, 5*points_within_1MHz, 5*points_within_1MHz) - lorentz = np.piecewise(x_filter, [(x_filter >= 0) * (x_filter < len(x_filter)*1/5), - (x_filter >= len(x_filter)*1/5) * (x_filter < len(x_filter)*2/5), - (x_filter >= len(x_filter)*2/5) * (x_filter < len(x_filter)*3/5), - (x_filter >= len(x_filter)*3/5) * (x_filter < len(x_filter)*4/5), - (x_filter >= len(x_filter)*4/5)], + x_filter = np.linspace(0, 5 * points_within_1MHz, 5 * points_within_1MHz) + lorentz = np.piecewise(x_filter, + [(x_filter >= 0) * (x_filter < len(x_filter) * 1 / 5), + (x_filter >= len(x_filter) * 1 / 5) * ( + x_filter < len(x_filter) * 2 / 5), + (x_filter >= len(x_filter) * 2 / 5) * ( + x_filter < len(x_filter) * 3 / 5), + (x_filter >= len(x_filter) * 3 / 5) * ( + x_filter < len(x_filter) * 4 / 5), + (x_filter >= len(x_filter) * 4 / 5)], [1, 0, 1, 0, 1]) # if the filter is smaller than 5 points a convolution does not make sense if len(lorentz) >= 5: data_convolved = filters.convolve1d(data_smooth_lorentz, - lorentz/lorentz.sum(), + lorentz / lorentz.sum(), mode='constant', cval=data_smooth_lorentz.max()) - x_axis_min = x_axis[data_convolved.argmin()]-2.15*1e6 + x_axis_min = x_axis[data_convolved.argmin()] - 2.15 * 1e6 else: - x_axis_min = x_axis[data_smooth_lorentz.argmin()]-2.15*1e6 + x_axis_min = x_axis[data_smooth_lorentz.argmin()] - 2.15 * 1e6 # level of the data, that means the offset is subtracted and the real data # are present @@ -921,23 +929,22 @@ def estimate_lorentziantriple_N14(self, x_axis, data, params): # That increases the accuracy of the calculated Integral. # integral of data corresponds to sqrt(2) * Amplitude * Sigma - smoothing_spline = 1 # must be 1<= smoothing_spline <= 5 + smoothing_spline = 1 # must be 1<= smoothing_spline <= 5 fit_function = InterpolatedUnivariateSpline(x_axis, data_level, k=smoothing_spline) integrated_area = fit_function.integral(x_axis[0], x_axis[-1]) - # sigma = abs(integrated_area / (minimum_level/np.pi)) # That is wrong, so commenting out: - sigma = abs(integrated_area /(np.pi * minimum_level))/3 + sigma = abs(integrated_area / (np.pi * minimum_level)) / 3 - amplitude = -1*abs(minimum_level) + amplitude = -1 * abs(minimum_level) # Since the total amplitude of the lorentzian is depending on sigma it makes # sense to vary sigma within an interval, which is smaller than the minimal # distance between two points. Then the fit algorithm will have a larger # range to determine the amplitude properly. That is the main issue with the # fit! - minimal_linewidth = (x_axis[1]-x_axis[0])/4 - maximal_linewidth = x_axis[-1]-x_axis[0] + minimal_linewidth = (x_axis[1] - x_axis[0]) / 4 + maximal_linewidth = x_axis[-1] - x_axis[0] # The linewidth of all the lorentzians are set to be the same! that is a # physical constraint for the N14 fitting. @@ -949,13 +956,13 @@ def estimate_lorentziantriple_N14(self, x_axis, data, params): params['l0_sigma'].set(value=sigma, min=minimal_linewidth, max=maximal_linewidth) params['l1_amplitude'].set(value=amplitude, max=-1e-6) - params['l1_center'].set(value=x_axis_min+hf_splitting, - expr='l0_center+{0}'.format(hf_splitting)) + params['l1_center'].set(value=x_axis_min + hf_splitting, + expr=f'l0_center+{hf_splitting}') params['l1_sigma'].set(value=sigma, min=minimal_linewidth, max=maximal_linewidth, expr='l0_sigma') params['l2_amplitude'].set(value=amplitude, max=-1e-6) - params['l2_center'].set(value=x_axis_min+hf_splitting*2, - expr='l0_center+{0}'.format(hf_splitting*2)) + params['l2_center'].set(value=x_axis_min + hf_splitting * 2, + expr=f'l0_center+{hf_splitting * 2}') params['l2_sigma'].set(value=sigma, min=minimal_linewidth, max=maximal_linewidth, expr='l0_sigma') params['offset'].set(value=offset) diff --git a/qudi_hira_analysis/_fitmethods/poissonianlikemethods.py b/qudi_hira_analysis/_fitmethods/poissonianlikemethods.py index 332a957..7d9fe48 100644 --- a/qudi_hira_analysis/_fitmethods/poissonianlikemethods.py +++ b/qudi_hira_analysis/_fitmethods/poissonianlikemethods.py @@ -1,5 +1,3 @@ -# -*- coding: utf-8 -*- - """ This file contains the Qudi fitting logic functions needed for poissinian-like-methods. @@ -24,14 +22,12 @@ from collections import OrderedDict import numpy as np -from lmfit import Parameters from lmfit.models import Model from scipy.interpolate import InterpolatedUnivariateSpline from scipy.ndimage import filters from scipy.signal import gaussian from scipy.special import gammaln, xlogy - ################################################################################ # # # Defining Poissonian models # @@ -49,10 +45,7 @@ def poisson(self, x, mu): Author: Travis Oliphant 2002-2011 with contributions from SciPy Developers 2004-2011 """ - if len(np.atleast_1d(x)) == 1: - check_val = x - else: - check_val = x[0] + check_val = x if len(np.atleast_1d(x)) == 1 else x[0] if check_val > 1e18: self.log.warning('The current value in the poissonian distribution ' @@ -94,6 +87,7 @@ def make_poissonian_model(self, prefix=None): lmfit.parameter.Parameter (without s) objects, keeping the information about the current value. """ + def poisson_function(x, mu): """ Function of a poisson distribution. @@ -108,7 +102,7 @@ def poisson_function(x, mu): if not isinstance(prefix, str) and prefix is not None: - self.log.error('The passed prefix <{0}> of type {1} is not a string and' + self.log.error('The passed prefix <{}> of type {} is not a string and' 'cannot be used as a prefix and will be ignored for now.' 'Correct that!'.format(prefix, type(prefix))) @@ -141,14 +135,17 @@ def make_poissonianmultiple_model(self, no_of_functions=1): multi_poisson_model, params = self.make_poissonian_model(prefix='p0_') for ii in range(1, no_of_functions): - multi_poisson_model += self.make_poissonian_model(prefix='p{0:d}_'.format(ii))[0] + multi_poisson_model += \ + self.make_poissonian_model(prefix=f'p{ii:d}_')[0] params = multi_poisson_model.make_params() return multi_poisson_model, params + def make_poissoniandouble_model(self): return self.make_poissonianmultiple_model(2) + ################################################################################ # # # Poissonian fits and their estimators # @@ -156,7 +153,8 @@ def make_poissoniandouble_model(self): ################################################################################ -def make_poissonian_fit(self, x_axis, data, estimator, units=None, add_params=None, **kwargs): +def make_poissonian_fit(self, x_axis, data, estimator, units=None, add_params=None, + **kwargs): """ Performe a poissonian fit on the provided data. @param numpy.array x_axis: 1D axis values @@ -193,15 +191,15 @@ def make_poissonian_fit(self, x_axis, data, estimator, units=None, add_params=No if units is None: units = ['arb. unit', 'arb. unit'] - result_str_dict = dict() # create result string for gui oder OrderedDict() + result_str_dict = {} # create result string for gui oder OrderedDict() result_str_dict['Amplitude'] = {'value': result.params['amplitude'].value, 'error': result.params['amplitude'].stderr, - 'unit': units[1]} # Amplitude + 'unit': units[1]} # Amplitude result_str_dict['Event rate'] = {'value': result.params['mu'].value, - 'error': result.params['mu'].stderr, - 'unit': units[0]} # event rate + 'error': result.params['mu'].stderr, + 'unit': units[0]} # event rate result.result_str_dict = result_str_dict @@ -227,8 +225,6 @@ def estimate_poissonian(self, x_axis, data, params): # a gaussian filter is appropriate due to the well approximation of poisson # distribution - # gaus = gaussian(10,10) - # data_smooth = filters.convolve1d(data, gaus/gaus.sum(), mode='mirror') data_smooth = self.gaussian_smoothing(data=data, filter_len=10, filter_sigma=10) @@ -240,7 +236,8 @@ def estimate_poissonian(self, x_axis, data, params): return error, params -def make_poissoniandouble_fit(self, x_axis, data, estimator, units=None, add_params=None, **kwargs): +def make_poissoniandouble_fit(self, x_axis, data, estimator, units=None, + add_params=None, **kwargs): """ Perform a double poissonian fit on the provided data. @param numpy.array x_axis: 1D axis values @@ -285,7 +282,7 @@ def make_poissoniandouble_fit(self, x_axis, data, estimator, units=None, add_par result_str_dict['Event rate 1'] = {'value': result.params['p0_mu'].value, 'error': result.params['p0_mu'].stderr, - 'unit': units[1]} + 'unit': units[1]} result_str_dict['Amplitude 2'] = {'value': result.params['p1_amplitude'].value, 'error': result.params['p1_amplitude'].stderr, @@ -293,7 +290,7 @@ def make_poissoniandouble_fit(self, x_axis, data, estimator, units=None, add_par result_str_dict['Event rate 2'] = {'value': result.params['p1_mu'].value, 'error': result.params['p1_mu'].stderr, - 'unit': units[1]} + 'unit': units[1]} result.result_str_dict = result_str_dict @@ -345,10 +342,7 @@ def estimate_poissoniandouble(self, x_axis, data, params, threshold_fraction=0.4 len_x = 10 interpol_factor = 1 else: - if len(x_axis) < 60: - interpol_factor = 4 - else: - interpol_factor = 2 + interpol_factor = 4 if len(x_axis) < 60 else 2 len_x = int(len(x_axis) / 10.) + 1 # Create the interpolation function, based on the data: @@ -376,11 +370,13 @@ def estimate_poissoniandouble(self, x_axis, data, params, threshold_fraction=0.4 # set the initial values for the fit: params['p0_mu'].set(value=x_axis_interpol[dip0_arg]) - amplitude0 = (data_smooth[dip0_arg] / self.poisson(x_axis_interpol[dip0_arg], x_axis_interpol[dip0_arg])) + amplitude0 = (data_smooth[dip0_arg] / self.poisson(x_axis_interpol[dip0_arg], + x_axis_interpol[dip0_arg])) params['p0_amplitude'].set(value=amplitude0, min=1e-15) params['p1_mu'].set(value=x_axis_interpol[dip1_arg]) - amplitude1 = (data_smooth[dip1_arg] / self.poisson(x_axis_interpol[dip1_arg], x_axis_interpol[dip1_arg])) + amplitude1 = (data_smooth[dip1_arg] / self.poisson(x_axis_interpol[dip1_arg], + x_axis_interpol[dip1_arg])) params['p1_amplitude'].set(value=amplitude1, min=1e-15) return error, params diff --git a/qudi_hira_analysis/_fitmethods/sinemethods.py b/qudi_hira_analysis/_fitmethods/sinemethods.py index 3043c9d..9448318 100644 --- a/qudi_hira_analysis/_fitmethods/sinemethods.py +++ b/qudi_hira_analysis/_fitmethods/sinemethods.py @@ -1,4 +1,3 @@ -# -*- coding: utf-8 -*- """ This file contains methods for sine fitting, these methods are imported by class FitLogic. @@ -42,17 +41,18 @@ def get_ft_windows(): """ win = {'none': {'func': np.ones, 'ampl_norm': 1.0}, - 'hamming': {'func': signal.hamming, 'ampl_norm': 1.0/0.54}, - 'hann': {'func': signal.hann, 'ampl_norm': 1.0/0.5}, - 'blackman': {'func': signal.blackman, 'ampl_norm': 1.0/0.42}, - 'triang': {'func': signal.triang, 'ampl_norm': 1.0/0.5}, - 'flattop': {'func': signal.flattop, 'ampl_norm': 1.0/0.2156}, - 'bartlett': {'func': signal.bartlett, 'ampl_norm': 1.0/0.5}, - 'parzen': {'func': signal.parzen, 'ampl_norm': 1.0/0.375}, - 'bohman': {'func': signal.bohman, 'ampl_norm': 1.0/0.4052847}, - 'blackmanharris': {'func': signal.blackmanharris, 'ampl_norm': 1.0/0.35875}, - 'nuttall': {'func': signal.nuttall, 'ampl_norm': 1.0/0.3635819}, - 'barthann': {'func': signal.barthann, 'ampl_norm': 1.0/0.5}} + 'hamming': {'func': signal.hamming, 'ampl_norm': 1.0 / 0.54}, + 'hann': {'func': signal.hann, 'ampl_norm': 1.0 / 0.5}, + 'blackman': {'func': signal.blackman, 'ampl_norm': 1.0 / 0.42}, + 'triang': {'func': signal.triang, 'ampl_norm': 1.0 / 0.5}, + 'flattop': {'func': signal.flattop, 'ampl_norm': 1.0 / 0.2156}, + 'bartlett': {'func': signal.bartlett, 'ampl_norm': 1.0 / 0.5}, + 'parzen': {'func': signal.parzen, 'ampl_norm': 1.0 / 0.375}, + 'bohman': {'func': signal.bohman, 'ampl_norm': 1.0 / 0.4052847}, + 'blackmanharris': {'func': signal.blackmanharris, + 'ampl_norm': 1.0 / 0.35875}, + 'nuttall': {'func': signal.nuttall, 'ampl_norm': 1.0 / 0.3635819}, + 'barthann': {'func': signal.barthann, 'ampl_norm': 1.0 / 0.5}} return win @@ -131,7 +131,7 @@ def compute_ft(x_val, y_val, zeropad_num=0, window='none', base_corr=True, psd=F ampl_norm_fact = avail_windows[window]['ampl_norm'] # zeropad for sinc interpolation: - zeropad_arr = np.zeros(len(corrected_y)*(zeropad_num+1)) + zeropad_arr = np.zeros(len(corrected_y) * (zeropad_num + 1)) zeropad_arr[:len(corrected_y)] = corrected_y # Get the amplitude values from the fourier transformed y values. @@ -145,13 +145,13 @@ def compute_ft(x_val, y_val, zeropad_num=0, window='none', base_corr=True, psd=F # The factor 2 accounts for the fact that just the half of the spectrum was # taken. The ampl_norm_fact is the normalization factor due to the applied # window function (the offset value in the window function): - fft_y = ((2/len(y_val)) * fft_y * ampl_norm_fact)**power_value + fft_y = ((2 / len(y_val)) * fft_y * ampl_norm_fact) ** power_value # Due to the sampling theorem you can only identify frequencies at half # of the sample rate, therefore the FT contains an almost symmetric # spectrum (the asymmetry results from aliasing effects). Therefore take # the half of the values for the display. - middle = int((len(zeropad_arr)+1)//2) + middle = int((len(zeropad_arr) + 1) // 2) # sample spacing of x_axis, if x is a time axis than it corresponds to a # timestep: @@ -210,10 +210,10 @@ def bare_sine_function(x, frequency, phase): model """ - return np.sin(2*np.pi*frequency*x+phase) + return np.sin(2 * np.pi * frequency * x + phase) if not isinstance(prefix, str) and prefix is not None: - self.log.error('The passed prefix <{0}> of type {1} is not a string and' + self.log.error('The passed prefix <{}> of type {} is not a string and' 'cannot be used as a prefix and will be ignored for now.' 'Correct that!'.format(prefix, type(prefix))) model = Model(bare_sine_function, independent_vars=['x']) @@ -224,6 +224,7 @@ def bare_sine_function(x, frequency, phase): return model, params + ############################### # Centred sine with no offset # ############################### @@ -242,11 +243,12 @@ def make_sinewithoutoffset_model(self, prefix=None): baresine_model, params = self.make_baresine_model(prefix=prefix) amplitude_model, params = self.make_amplitude_model(prefix=prefix) - sine_model = amplitude_model*baresine_model + sine_model = amplitude_model * baresine_model params = sine_model.make_params() return sine_model, params + ############################ # General sine with offset # ############################ @@ -271,6 +273,7 @@ def make_sine_model(self, prefix=None): return sine_offset_model, params + ############################################### # Sinus with exponential decay but not offset # ############################################### @@ -288,13 +291,15 @@ def make_sineexpdecaywithoutoffset_model(self, prefix=None): """ sine_model, params = self.make_sinewithoutoffset_model(prefix=prefix) - bareexponentialdecay_model, params = self.make_bareexponentialdecay_model(prefix=prefix) + bareexponentialdecay_model, params = self.make_bareexponentialdecay_model( + prefix=prefix) - sine_exp_decay_model = sine_model*bareexponentialdecay_model + sine_exp_decay_model = sine_model * bareexponentialdecay_model params = sine_exp_decay_model.make_params() return sine_exp_decay_model, params + ################################################### # Sinus with exponential decay and offset fitting # ################################################### @@ -310,7 +315,8 @@ def make_sineexponentialdecay_model(self, prefix=None): the method make_baresine_model. """ - sine_exp_decay_model, params = self.make_sineexpdecaywithoutoffset_model(prefix=prefix) + sine_exp_decay_model, params = self.make_sineexpdecaywithoutoffset_model( + prefix=prefix) constant_model, params = self.make_constant_model(prefix=prefix) sine_exp_decay_offset_model = sine_exp_decay_model + constant_model @@ -318,6 +324,7 @@ def make_sineexponentialdecay_model(self, prefix=None): return sine_exp_decay_offset_model, params + ################################################### # Sinus with stretched exponential decay fitting # ################################################### @@ -334,7 +341,8 @@ def make_sinestretchedexponentialdecay_model(self, prefix=None): """ sine_model, params = self.make_sinewithoutoffset_model(prefix=prefix) - bare_stretched_exp_decay_model, params = self.make_barestretchedexponentialdecay_model(prefix=prefix) + bare_stretched_exp_decay_model, params = self.make_barestretchedexponentialdecay_model( + prefix=prefix) constant_model, params = self.make_constant_model(prefix=prefix) model = sine_model * bare_stretched_exp_decay_model + constant_model @@ -342,6 +350,7 @@ def make_sinestretchedexponentialdecay_model(self, prefix=None): return model, params + ########################################### # Sum of two individual Sinus with offset # ########################################### @@ -358,13 +367,10 @@ def make_sinedouble_model(self, prefix=None): the method make_baresine_model. """ - if prefix is None: - add_text = '' - else: - add_text = prefix + add_text = '' if prefix is None else prefix - sine_model1, params = self.make_sinewithoutoffset_model(prefix='s1_'+add_text) - sine_model2, params = self.make_sinewithoutoffset_model(prefix='s2_'+add_text) + sine_model1, params = self.make_sinewithoutoffset_model(prefix='s1_' + add_text) + sine_model2, params = self.make_sinewithoutoffset_model(prefix='s2_' + add_text) constant_model, params = self.make_constant_model(prefix=prefix) @@ -373,6 +379,7 @@ def make_sinedouble_model(self, prefix=None): return two_sine_offset, params + ################################################################################ # Sum of two individual Sinus with offset and single exponential decay # ################################################################################ @@ -389,22 +396,21 @@ def make_sinedoublewithexpdecay_model(self, prefix=None): the method make_baresine_model. """ - if prefix is None: - add_text = '' - else: - add_text = prefix + add_text = '' if prefix is None else prefix - sine_model1, params = self.make_sinewithoutoffset_model(prefix='s1_'+add_text) - sine_model2, params = self.make_sinewithoutoffset_model(prefix='s2_'+add_text) + sine_model1, params = self.make_sinewithoutoffset_model(prefix='s1_' + add_text) + sine_model2, params = self.make_sinewithoutoffset_model(prefix='s2_' + add_text) bare_exp_decay_model, params = self.make_bareexponentialdecay_model(prefix=prefix) constant_model, params = self.make_constant_model(prefix=prefix) - two_sine_exp_decay_offset = (sine_model1 + sine_model2)*bare_exp_decay_model + constant_model + two_sine_exp_decay_offset = ( + sine_model1 + sine_model2) * bare_exp_decay_model + constant_model params = two_sine_exp_decay_offset.make_params() return two_sine_exp_decay_offset, params + ############################################################### # Sum of two individual Sinus exponential decays (and offset) # ############################################################### @@ -420,13 +426,12 @@ def make_sinedoublewithtwoexpdecay_model(self, prefix=None): @return tuple: (object model, object params), for more description see in the method make_baresine_model. """ - if prefix is None: - add_text = '' - else: - add_text = prefix + add_text = '' if prefix is None else prefix - sine_exp_decay_model1, params = self.make_sineexpdecaywithoutoffset_model(prefix='e1_'+add_text) - sine_exp_decay_model2, params = self.make_sineexpdecaywithoutoffset_model(prefix='e2_'+add_text) + sine_exp_decay_model1, params = self.make_sineexpdecaywithoutoffset_model( + prefix='e1_' + add_text) + sine_exp_decay_model2, params = self.make_sineexpdecaywithoutoffset_model( + prefix='e2_' + add_text) constant_model, params = self.make_constant_model(prefix=prefix) @@ -435,6 +440,7 @@ def make_sinedoublewithtwoexpdecay_model(self, prefix=None): return sinedoublewithtwoexpdecay, params + ############################################# # Sum of three individual Sinus with offset # ############################################# @@ -451,14 +457,11 @@ def make_sinetriple_model(self, prefix=None): the method make_baresine_model. """ - if prefix is None: - add_text = '' - else: - add_text = prefix + add_text = '' if prefix is None else prefix - sine_model1, params = self.make_sinewithoutoffset_model(prefix='s1_'+add_text) - sine_model2, params = self.make_sinewithoutoffset_model(prefix='s2_'+add_text) - sine_model3, params = self.make_sinewithoutoffset_model(prefix='s3_'+add_text) + sine_model1, params = self.make_sinewithoutoffset_model(prefix='s1_' + add_text) + sine_model2, params = self.make_sinewithoutoffset_model(prefix='s2_' + add_text) + sine_model3, params = self.make_sinewithoutoffset_model(prefix='s3_' + add_text) constant_model, params = self.make_constant_model(prefix=prefix) @@ -467,6 +470,7 @@ def make_sinetriple_model(self, prefix=None): return three_sine_offset, params + ########################################################################## # Sum of three individual Sinus with offset and single exponential decay # ########################################################################## @@ -483,24 +487,22 @@ def make_sinetriplewithexpdecay_model(self, prefix=None): the method make_baresine_model. """ - if prefix is None: - add_text = '' - else: - add_text = prefix + add_text = '' if prefix is None else prefix - sine_model1, params = self.make_sinewithoutoffset_model(prefix='s1_'+add_text) - sine_model2, params = self.make_sinewithoutoffset_model(prefix='s2_'+add_text) - sine_model3, params = self.make_sinewithoutoffset_model(prefix='s3_'+add_text) + sine_model1, params = self.make_sinewithoutoffset_model(prefix='s1_' + add_text) + sine_model2, params = self.make_sinewithoutoffset_model(prefix='s2_' + add_text) + sine_model3, params = self.make_sinewithoutoffset_model(prefix='s3_' + add_text) bare_exp_decay_model, params = self.make_bareexponentialdecay_model(prefix=prefix) constant_model, params = self.make_constant_model(prefix=prefix) three_sine_exp_decay_offset = ( - sine_model1 + sine_model2 + sine_model3) * bare_exp_decay_model + constant_model + sine_model1 + sine_model2 + sine_model3) * bare_exp_decay_model + constant_model params = three_sine_exp_decay_offset.make_params() return three_sine_exp_decay_offset, params + ######################################################################### # Sum of three individual Sinus with offset and three exponential decay # ######################################################################### @@ -517,14 +519,14 @@ def make_sinetriplewiththreeexpdecay_model(self, prefix=None): the method make_baresine_model. """ - if prefix is None: - add_text = '' - else: - add_text = prefix + add_text = '' if prefix is None else prefix - sine_exp_decay_model1, params = self.make_sineexpdecaywithoutoffset_model(prefix='e1_'+add_text) - sine_exp_decay_model2, params = self.make_sineexpdecaywithoutoffset_model(prefix='e2_'+add_text) - sine_exp_decay_model3, params = self.make_sineexpdecaywithoutoffset_model(prefix='e3_'+add_text) + sine_exp_decay_model1, params = self.make_sineexpdecaywithoutoffset_model( + prefix='e1_' + add_text) + sine_exp_decay_model2, params = self.make_sineexpdecaywithoutoffset_model( + prefix='e2_' + add_text) + sine_exp_decay_model3, params = self.make_sineexpdecaywithoutoffset_model( + prefix='e3_' + add_text) constant_model, params = self.make_constant_model(prefix=prefix) @@ -533,6 +535,7 @@ def make_sinetriplewiththreeexpdecay_model(self, prefix=None): return three_sine_exp_decay_offset, params + ################################################################################ # # # General estimators used later # @@ -566,14 +569,14 @@ def estimate_baresine(self, x_axis, data, params): # appearing peaks. dft_x, dft_y = compute_ft(x_axis, data, zeropad_num=1) - stepsize = x_axis[1]-x_axis[0] # for frequency axis + stepsize = x_axis[1] - x_axis[0] # for frequency axis frequency_max = np.abs(dft_x[np.log(dft_y).argmax()]) # find minimal distance to the next meas point in the corresponding time value> min_x_diff = np.ediff1d(x_axis).min() # How many points are used to sample the estimated frequency with min_x_diff: - iter_steps = int(1/(frequency_max*min_x_diff)) + iter_steps = int(1 / (frequency_max * min_x_diff)) if iter_steps < 1: iter_steps = 1 @@ -585,15 +588,16 @@ def estimate_baresine(self, x_axis, data, params): # trace. for iter_s in range(iter_steps): - func_val = np.sin(2*np.pi*frequency_max*x_axis + iter_s/iter_steps *2*np.pi) + func_val = np.sin( + 2 * np.pi * frequency_max * x_axis + iter_s / iter_steps * 2 * np.pi) sum_res[iter_s] = np.abs(data - func_val).sum() # The minimum indicates where the sine function was fittng the worst, # therefore subtract pi. This will also ensure that the estimated phase will # be in the interval [-pi,pi]. - phase = sum_res.argmax()/iter_steps *2*np.pi - np.pi + phase = sum_res.argmax() / iter_steps * 2 * np.pi - np.pi - params['frequency'].set(value=frequency_max, min=0.0, max=1/stepsize*3) + params['frequency'].set(value=frequency_max, min=0.0, max=1 / stepsize * 3) params['phase'].set(value=phase, min=-np.pi, max=np.pi) return error, params @@ -652,7 +656,7 @@ def estimate_sinewithoutoffset(self, x_axis, data, params): # if at least two identical values are in the array, then the difference is of course zero, # catch that case. - for tries, diff_array_step in enumerate(diff_array): + for _tries, _diff_array_step in enumerate(diff_array): if np.isclose(min_x_diff, 0.0, atol=1e-12): index = np.argmin(diff_array) diff_array = np.delete(diff_array, index) @@ -670,7 +674,7 @@ def estimate_sinewithoutoffset(self, x_axis, data, params): break # How many points are used to sample the estimated frequency with min_x_diff: - iter_steps = int(1/(frequency_max*min_x_diff)) + iter_steps = int(1 / (frequency_max * min_x_diff)) if iter_steps < 1: iter_steps = 1 @@ -682,17 +686,18 @@ def estimate_sinewithoutoffset(self, x_axis, data, params): # trace. for iter_s in range(iter_steps): - func_val = ampl_val * np.sin(2*np.pi*frequency_max*x_axis + iter_s/iter_steps*2*np.pi) + func_val = ampl_val * np.sin( + 2 * np.pi * frequency_max * x_axis + iter_s / iter_steps * 2 * np.pi) sum_res[iter_s] = np.abs(data - func_val).sum() # The minimum indicates where the sine function was fitting the worst, # therefore subtract pi. This will also ensure that the estimated phase will # be in the interval [-pi,pi]. - phase = sum_res.argmax()/iter_steps *2*np.pi - np.pi + phase = sum_res.argmax() / iter_steps * 2 * np.pi - np.pi # values and bounds of initial parameters params['amplitude'].set(value=ampl_val) - params['frequency'].set(value=frequency_max, min=0.0, max=1/stepsize*3) + params['frequency'].set(value=frequency_max, min=0.0, max=1 / stepsize * 3) params['phase'].set(value=phase, min=-np.pi, max=np.pi) return error, params @@ -736,16 +741,17 @@ def make_sine_fit(self, x_axis, data, estimator, units=None, add_params=None, ** except: result = sine.fit(data, x=x_axis, params=params, **kwargs) self.log.error('The sine fit did not work.\n' - 'Error message: {0}\n'.format(result.message)) + f'Error message: {result.message}\n') if units is None: units = ['arb. unit', 'arb. unit'] - result_str_dict = dict() + result_str_dict = {} period = 1 / result.params['frequency'].value try: - period_err = result.params['frequency'].stderr / (result.params['frequency'])**2 + period_err = result.params['frequency'].stderr / ( + result.params['frequency']) ** 2 except ZeroDivisionError: period_err = np.inf @@ -755,10 +761,11 @@ def make_sine_fit(self, x_axis, data, estimator, units=None, add_params=None, ** result_str_dict['Frequency'] = {'value': result.params['frequency'].value, 'error': result.params['frequency'].stderr, - 'unit': 'Hz' if units[0] == 's' else '1/' + units[0]} + 'unit': 'Hz' if units[0] == 's' else '1/' + units[ + 0]} - result_str_dict['Phase'] = {'value': 180/np.pi*result.params['phase'].value, - 'error': 180/np.pi*result.params['phase'].stderr, + result_str_dict['Phase'] = {'value': 180 / np.pi * result.params['phase'].value, + 'error': 180 / np.pi * result.params['phase'].stderr, 'unit': 'deg'} result_str_dict['Offset'] = {'value': result.params['offset'].value, @@ -766,18 +773,28 @@ def make_sine_fit(self, x_axis, data, estimator, units=None, add_params=None, ** 'unit': units[1]} result_str_dict['Amplitude'] = {'value': result.params['amplitude'].value, - 'error': result.params['amplitude'].stderr, + 'error': result.params['amplitude'].stderr, 'unit': units[1]} - result_str_dict['Contrast'] = {'value': ((2*result.params['amplitude'].value) / - (result.params['offset'].value+result.params['amplitude'].value)*100), - 'error': (np.abs((2*result.params['amplitude'].value) / - (result.params['offset'].value+result.params['amplitude'].value)**2 * - result.params['offset'].stderr) + - np.abs((2/(result.params['offset'].value + - result.params['amplitude'].value) + (2*result.params['amplitude'].value) / - (result.params['offset'].value + result.params['amplitude'].value)**2) * - result.params['amplitude'].stderr))*100, + result_str_dict['Contrast'] = {'value': ((2 * result.params['amplitude'].value) / + (result.params['offset'].value + + result.params['amplitude'].value) * 100), + 'error': (np.abs( + (2 * result.params['amplitude'].value) / + (result.params['offset'].value + result.params[ + 'amplitude'].value) ** 2 * + result.params['offset'].stderr) + + np.abs( + (2 / (result.params['offset'].value + + result.params[ + 'amplitude'].value) + ( + 2 * result.params[ + 'amplitude'].value) / + (result.params['offset'].value + + result.params[ + 'amplitude'].value) ** 2) * + result.params[ + 'amplitude'].stderr)) * 100, 'unit': '%'} result.result_str_dict = result_str_dict @@ -807,18 +824,21 @@ def estimate_sine(self, x_axis, data, params): # level data data_level = data - offset - error, params = self.estimate_sinewithoutoffset(x_axis=x_axis, data=data_level, params=params) + error, params = self.estimate_sinewithoutoffset(x_axis=x_axis, data=data_level, + params=params) params['offset'].set(value=offset) return error, params + ########################## # Sine exponential decay # ########################## -def make_sineexponentialdecay_fit(self, x_axis, data, estimator, units=None, add_params=None, **kwargs): +def make_sineexponentialdecay_fit(self, x_axis, data, estimator, units=None, + add_params=None, **kwargs): """ Perform a sine exponential decay fit on the provided data. @param numpy.array x_axis: 1D axis values @@ -846,16 +866,17 @@ def make_sineexponentialdecay_fit(self, x_axis, data, estimator, units=None, add result = sine_exp_decay_offset.fit(data, x=x_axis, params=params, **kwargs) self.log.error('The sineexponentialdecayoffset fit did not work.\n' - 'Error message: {0}'.format(result.message)) + f'Error message: {result.message}') if units is None: units = ['arb. unit', 'arb. unit'] - result_str_dict = dict() + result_str_dict = {} - period = 1/result.params['frequency'].value + period = 1 / result.params['frequency'].value try: - period_err = result.params['frequency'].stderr / (result.params['frequency'])**2 + period_err = result.params['frequency'].stderr / ( + result.params['frequency']) ** 2 except ZeroDivisionError: period_err = np.inf @@ -865,25 +886,36 @@ def make_sineexponentialdecay_fit(self, x_axis, data, estimator, units=None, add result_str_dict['Frequency'] = {'value': result.params['frequency'].value, 'error': result.params['frequency'].stderr, - 'unit': 'Hz' if units[0] == 's' else '1/'+units[0]} + 'unit': 'Hz' if units[0] == 's' else '1/' + units[ + 0]} result_str_dict['Amplitude'] = {'value': result.params['amplitude'].value, 'error': result.params['amplitude'].stderr, 'unit': units[1]} - result_str_dict['Contrast'] = {'value': ((2*result.params['amplitude'].value) / - (result.params['offset'].value+result.params['amplitude'].value)*100), - 'error': (np.abs((2*result.params['amplitude'].value) / - (result.params['offset'].value+result.params['amplitude'].value)**2 * - result.params['offset'].stderr) + - np.abs((2/(result.params['offset'].value + - result.params['amplitude'].value) + (2*result.params['amplitude'].value) / - (result.params['offset'].value + result.params['amplitude'].value)**2) * - result.params['amplitude'].stderr))*100, + result_str_dict['Contrast'] = {'value': ((2 * result.params['amplitude'].value) / + (result.params['offset'].value + + result.params['amplitude'].value) * 100), + 'error': (np.abs( + (2 * result.params['amplitude'].value) / + (result.params['offset'].value + result.params[ + 'amplitude'].value) ** 2 * + result.params['offset'].stderr) + + np.abs( + (2 / (result.params['offset'].value + + result.params[ + 'amplitude'].value) + ( + 2 * result.params[ + 'amplitude'].value) / + (result.params['offset'].value + + result.params[ + 'amplitude'].value) ** 2) * + result.params[ + 'amplitude'].stderr)) * 100, 'unit': '%'} - result_str_dict['Phase'] = {'value': 180/np.pi*result.params['phase'].value, - 'error': 180/np.pi*result.params['phase'].stderr, + result_str_dict['Phase'] = {'value': 180 / np.pi * result.params['phase'].value, + 'error': 180 / np.pi * result.params['phase'].stderr, 'unit': 'deg'} result_str_dict['Offset'] = {'value': result.params['offset'].value, @@ -941,21 +973,21 @@ def estimate_sineexponentialdecay(self, x_axis, data, params=None): # remove noise a = np.std(dft_y) - for i in range(0, len(dft_x)): + for i in range(len(dft_x)): if dft_y[i] <= a: dft_y[i] = 0 # calculating the width of the FT peak for the estimation of lifetime s = 0 - for i in range(0, len(dft_x)): - s += dft_y[i]*abs(dft_x[1]-dft_x[0])/max(dft_y) - lifetime_val = 0.5/s + for i in range(len(dft_x)): + s += dft_y[i] * abs(dft_x[1] - dft_x[0]) / max(dft_y) + lifetime_val = 0.5 / s # find minimal distance to the next meas point in the corresponding x value min_x_diff = np.ediff1d(x_axis).min() # How many points are used to sample the estimated frequency with min_x_diff: - iter_steps = int(1/(frequency_max*min_x_diff)) + iter_steps = int(1 / (frequency_max * min_x_diff)) if iter_steps < 1: iter_steps = 1 @@ -967,34 +999,38 @@ def estimate_sineexponentialdecay(self, x_axis, data, params=None): # trace. for iter_s in range(iter_steps): - func_val = ampl_val * np.sin(2*np.pi*frequency_max*x_axis + iter_s/iter_steps *2*np.pi) + func_val = ampl_val * np.sin( + 2 * np.pi * frequency_max * x_axis + iter_s / iter_steps * 2 * np.pi) sum_res[iter_s] = np.abs(data_level - func_val).sum() # The minimum indicates where the sine function was fittng the worst, # therefore subtract pi. This will also ensure that the estimated phase will # be in the interval [-pi,pi]. - phase = (sum_res.argmax()/iter_steps *2*np.pi - np.pi)%(2*np.pi) + phase = (sum_res.argmax() / iter_steps * 2 * np.pi - np.pi) % (2 * np.pi) # values and bounds of initial parameters params['frequency'].set(value=frequency_max, - min=min(0.1 / (x_axis[-1]-x_axis[0]), dft_x[3]), - max=min(0.5 / stepsize, dft_x.max()-abs(dft_x[2]-dft_x[0]))) - params['phase'].set(value=phase, min=-2*np.pi, max=2*np.pi) + min=min(0.1 / (x_axis[-1] - x_axis[0]), dft_x[3]), + max=min(0.5 / stepsize, + dft_x.max() - abs(dft_x[2] - dft_x[0]))) + params['phase'].set(value=phase, min=-2 * np.pi, max=2 * np.pi) params['amplitude'].set(value=ampl_val, min=0) params['offset'].set(value=offset) params['lifetime'].set(value=lifetime_val, - min=2*(x_axis[1]-x_axis[0]), - max=1/(abs(dft_x[1]-dft_x[0])*0.5)) + min=2 * (x_axis[1] - x_axis[0]), + max=1 / (abs(dft_x[1] - dft_x[0]) * 0.5)) return error, params + ################################################### # Sinus with stretched exponential decay fitting # ################################################### -def make_sinestretchedexponentialdecay_fit(self, x_axis, data, estimator, units=None, add_params=None, **kwargs): +def make_sinestretchedexponentialdecay_fit(self, x_axis, data, estimator, units=None, + add_params=None, **kwargs): """ Perform a sine stretched exponential decay fit on the provided data. @param numpy.array x_axis: 1D axis values @@ -1021,16 +1057,17 @@ def make_sinestretchedexponentialdecay_fit(self, x_axis, data, estimator, units= except: result = sine_stretched_exp_decay.fit(data, x=x_axis, params=params, **kwargs) self.log.error('The sineexponentialdecay fit did not work.\n' - 'Error message: {0}'.format(result.message)) + f'Error message: {result.message}') if units is None: units = ['arb. unit', 'arb. unit'] - result_str_dict = dict() + result_str_dict = {} period = 1 / result.params['frequency'].value try: - period_err = result.params['frequency'].stderr / (result.params['frequency'])**2 + period_err = result.params['frequency'].stderr / ( + result.params['frequency']) ** 2 except ZeroDivisionError: period_err = np.inf @@ -1040,25 +1077,36 @@ def make_sinestretchedexponentialdecay_fit(self, x_axis, data, estimator, units= result_str_dict['Frequency'] = {'value': result.params['frequency'].value, 'error': result.params['frequency'].stderr, - 'unit': 'Hz' if units[0] == 's' else '1/' + units[0]} + 'unit': 'Hz' if units[0] == 's' else '1/' + units[ + 0]} result_str_dict['Amplitude'] = {'value': result.params['amplitude'].value, 'error': result.params['amplitude'].stderr, 'unit': units[1]} - result_str_dict['Contrast'] = {'value': ((2*result.params['amplitude'].value) / - (result.params['offset'].value+result.params['amplitude'].value)*100), - 'error': (np.abs((2*result.params['amplitude'].value) / - (result.params['offset'].value+result.params['amplitude'].value)**2 * - result.params['offset'].stderr) + - np.abs((2/(result.params['offset'].value + - result.params['amplitude'].value) + (2*result.params['amplitude'].value) / - (result.params['offset'].value + result.params['amplitude'].value)**2) * - result.params['amplitude'].stderr))*100, + result_str_dict['Contrast'] = {'value': ((2 * result.params['amplitude'].value) / + (result.params['offset'].value + + result.params['amplitude'].value) * 100), + 'error': (np.abs( + (2 * result.params['amplitude'].value) / + (result.params['offset'].value + result.params[ + 'amplitude'].value) ** 2 * + result.params['offset'].stderr) + + np.abs( + (2 / (result.params['offset'].value + + result.params[ + 'amplitude'].value) + ( + 2 * result.params[ + 'amplitude'].value) / + (result.params['offset'].value + + result.params[ + 'amplitude'].value) ** 2) * + result.params[ + 'amplitude'].stderr)) * 100, 'unit': '%'} - result_str_dict['Phase'] = {'value': 180/np.pi*result.params['phase'].value, - 'error': 180/np.pi*result.params['phase'].stderr, + result_str_dict['Phase'] = {'value': 180 / np.pi * result.params['phase'].value, + 'error': 180 / np.pi * result.params['phase'].stderr, 'unit': 'deg'} result_str_dict['Offset'] = {'value': result.params['offset'].value, @@ -1092,18 +1140,20 @@ def estimate_sinestretchedexponentialdecay(self, x_axis, data, params): """ error, params = self.estimate_sineexponentialdecay(x_axis, data, params) - #TODO: estimate the exponent cleaverly! For now, set the initial value to 2 + # TODO: estimate the exponent cleaverly! For now, set the initial value to 2 # since the usual values for our cases are between 1 and 3. params['beta'].set(value=2, min=0.0, max=10) return error, params + ########################################### # Sum of two individual Sinus with offset # ########################################### -def make_sinedouble_fit(self, x_axis, data, estimator, units=None, add_params=None, **kwargs): +def make_sinedouble_fit(self, x_axis, data, estimator, units=None, add_params=None, + **kwargs): """ Perform a two sine with offset fit on the provided data. @param numpy.array x_axis: 1D axis values @@ -1129,23 +1179,25 @@ def make_sinedouble_fit(self, x_axis, data, estimator, units=None, add_params=No result = two_sine_offset.fit(data, x=x_axis, params=params, **kwargs) except: self.log.warning('The twosineexpdecayoffset fit did not work. ' - 'Error message: {}'.format(str(result.message))) + f'Error message: {result.message!s}') result = two_sine_offset.fit(data, x=x_axis, params=params, **kwargs) if units is None: units = ['arb. unit', 'arb. unit'] - result_str_dict = dict() # create result string for gui or OrderedDict() + result_str_dict = {} # create result string for gui or OrderedDict() period1 = 1 / result.params['s1_frequency'].value try: - period1_err = result.params['s1_frequency'].stderr / (result.params['s1_frequency'])**2 + period1_err = result.params['s1_frequency'].stderr / ( + result.params['s1_frequency']) ** 2 except ZeroDivisionError: period1_err = np.inf period2 = 1 / result.params['s2_frequency'].value try: - period2_err = result.params['s2_frequency'].stderr / (result.params['s2_frequency'])**2 + period2_err = result.params['s2_frequency'].stderr / ( + result.params['s2_frequency']) ** 2 except ZeroDivisionError: period2_err = np.inf @@ -1159,11 +1211,13 @@ def make_sinedouble_fit(self, x_axis, data, estimator, units=None, add_params=No result_str_dict['Frequency 1'] = {'value': result.params['s1_frequency'].value, 'error': result.params['s1_frequency'].stderr, - 'unit': 'Hz' if units[0] == 's' else '1/' + units[0]} + 'unit': 'Hz' if units[0] == 's' else '1/' + units[ + 0]} result_str_dict['Frequency 2'] = {'value': result.params['s2_frequency'].value, 'error': result.params['s2_frequency'].stderr, - 'unit': 'Hz' if units[0] == 's' else '1/' + units[0]} + 'unit': 'Hz' if units[0] == 's' else '1/' + units[ + 0]} result_str_dict['Amplitude 1'] = {'value': result.params['s1_amplitude'].value, 'error': result.params['s1_amplitude'].stderr, @@ -1174,36 +1228,58 @@ def make_sinedouble_fit(self, x_axis, data, estimator, units=None, add_params=No 'unit': units[1]} amp_string = 's1_amplitude' - result_str_dict['Contrast 1'] = {'value': ((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)*100), - 'error': (np.abs((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)**2 * - result.params['offset'].stderr) + - np.abs((2/(result.params['offset'].value + - result.params[amp_string].value) + (2*result.params[amp_string].value) / - (result.params['offset'].value + result.params[amp_string].value)**2) * - result.params[amp_string].stderr))*100, - 'unit': '%'} + result_str_dict['Contrast 1'] = {'value': ((2 * result.params[amp_string].value) / + (result.params['offset'].value + + result.params[amp_string].value) * 100), + 'error': (np.abs( + (2 * result.params[amp_string].value) / + (result.params['offset'].value + result.params[ + amp_string].value) ** 2 * + result.params['offset'].stderr) + + np.abs( + (2 / (result.params['offset'].value + + result.params[ + amp_string].value) + ( + 2 * result.params[ + amp_string].value) / + (result.params['offset'].value + + result.params[ + amp_string].value) ** 2) * + result.params[ + amp_string].stderr)) * 100, + 'unit': '%'} amp_string = 's2_amplitude' - result_str_dict['Contrast 2'] = {'value': ((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)*100), - 'error': (np.abs((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)**2 * - result.params['offset'].stderr) + - np.abs((2/(result.params['offset'].value + - result.params[amp_string].value) + (2*result.params[amp_string].value) / - (result.params['offset'].value + result.params[amp_string].value)**2) * - result.params[amp_string].stderr))*100, - 'unit': '%'} + result_str_dict['Contrast 2'] = {'value': ((2 * result.params[amp_string].value) / + (result.params['offset'].value + + result.params[amp_string].value) * 100), + 'error': (np.abs( + (2 * result.params[amp_string].value) / + (result.params['offset'].value + result.params[ + amp_string].value) ** 2 * + result.params['offset'].stderr) + + np.abs( + (2 / (result.params['offset'].value + + result.params[ + amp_string].value) + ( + 2 * result.params[ + amp_string].value) / + (result.params['offset'].value + + result.params[ + amp_string].value) ** 2) * + result.params[ + amp_string].stderr)) * 100, + 'unit': '%'} - result_str_dict['Phase 1'] = {'value': 180/np.pi*result.params['s1_phase'].value, - 'error': 180/np.pi*result.params['s1_phase'].stderr, - 'unit': 'deg'} + result_str_dict['Phase 1'] = { + 'value': 180 / np.pi * result.params['s1_phase'].value, + 'error': 180 / np.pi * result.params['s1_phase'].stderr, + 'unit': 'deg'} - result_str_dict['Phase 2'] = {'value': 180/np.pi*result.params['s2_phase'].value, - 'error': 180/np.pi*result.params['s2_phase'].stderr, - 'unit': 'deg'} + result_str_dict['Phase 2'] = { + 'value': 180 / np.pi * result.params['s2_phase'].value, + 'error': 180 / np.pi * result.params['s2_phase'].stderr, + 'unit': 'deg'} result_str_dict['Offset'] = {'value': result.params['offset'].value, 'error': result.params['offset'].stderr, @@ -1237,7 +1313,8 @@ def estimate_sinedouble(self, x_axis, data, params): result1 = self.make_sine_fit(x_axis=x_axis, data=data, estimator=self.estimate_sine) data_sub = data - result1.best_fit - result2 = self.make_sine_fit(x_axis=x_axis, data=data_sub, estimator=self.estimate_sine) + result2 = self.make_sine_fit(x_axis=x_axis, data=data_sub, + estimator=self.estimate_sine) # Fill the parameter dict: params['s1_amplitude'].set(value=result1.params['amplitude'].value) @@ -1252,12 +1329,14 @@ def estimate_sinedouble(self, x_axis, data, params): return error, params + ################################################################################ # Sum of two individual Sinus with offset and single exponential decay # ################################################################################ -def make_sinedoublewithexpdecay_fit(self, x_axis, data, estimator, units=None, add_params=None, **kwargs): +def make_sinedoublewithexpdecay_fit(self, x_axis, data, estimator, units=None, + add_params=None, **kwargs): """ Perform a two sine with one exponential decay offset fit on the provided data. @@ -1284,23 +1363,25 @@ def make_sinedoublewithexpdecay_fit(self, x_axis, data, estimator, units=None, a result = two_sine_exp_decay_offset.fit(data, x=x_axis, params=params, **kwargs) except: self.log.warning('The sinedoublewithexpdecay fit did not work. ' - 'Error message: {}'.format(str(result.message))) + f'Error message: {result.message!s}') result = two_sine_exp_decay_offset.fit(data, x=x_axis, params=params, **kwargs) if units is None: units = ['arb. unit', 'arb. unit'] - result_str_dict = dict() # create result string for gui or OrderedDict() + result_str_dict = {} # create result string for gui or OrderedDict() period1 = 1 / result.params['s1_frequency'].value try: - period1_err = result.params['s1_frequency'].stderr / (result.params['s1_frequency']) ** 2 + period1_err = result.params['s1_frequency'].stderr / ( + result.params['s1_frequency']) ** 2 except ZeroDivisionError: period1_err = np.inf period2 = 1 / result.params['s2_frequency'].value try: - period2_err = result.params['s2_frequency'].stderr / (result.params['s2_frequency']) ** 2 + period2_err = result.params['s2_frequency'].stderr / ( + result.params['s2_frequency']) ** 2 except ZeroDivisionError: period2_err = np.inf @@ -1314,11 +1395,13 @@ def make_sinedoublewithexpdecay_fit(self, x_axis, data, estimator, units=None, a result_str_dict['Frequency 1'] = {'value': result.params['s1_frequency'].value, 'error': result.params['s1_frequency'].stderr, - 'unit': 'Hz' if units[0] == 's' else '1/' + units[0]} + 'unit': 'Hz' if units[0] == 's' else '1/' + units[ + 0]} result_str_dict['Frequency 2'] = {'value': result.params['s2_frequency'].value, 'error': result.params['s2_frequency'].stderr, - 'unit': 'Hz' if units[0] == 's' else '1/' + units[0]} + 'unit': 'Hz' if units[0] == 's' else '1/' + units[ + 0]} result_str_dict['Amplitude 1'] = {'value': result.params['s1_amplitude'].value, 'error': result.params['s1_amplitude'].stderr, @@ -1329,36 +1412,58 @@ def make_sinedoublewithexpdecay_fit(self, x_axis, data, estimator, units=None, a 'unit': units[1]} amp_string = 's1_amplitude' - result_str_dict['Contrast 1'] = {'value': ((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)*100), - 'error': (np.abs((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)**2 * - result.params['offset'].stderr) + - np.abs((2/(result.params['offset'].value + - result.params[amp_string].value) + (2*result.params[amp_string].value) / - (result.params['offset'].value + result.params[amp_string].value)**2) * - result.params[amp_string].stderr))*100, + result_str_dict['Contrast 1'] = {'value': ((2 * result.params[amp_string].value) / + (result.params['offset'].value + + result.params[amp_string].value) * 100), + 'error': (np.abs( + (2 * result.params[amp_string].value) / + (result.params['offset'].value + result.params[ + amp_string].value) ** 2 * + result.params['offset'].stderr) + + np.abs( + (2 / (result.params['offset'].value + + result.params[ + amp_string].value) + ( + 2 * result.params[ + amp_string].value) / + (result.params['offset'].value + + result.params[ + amp_string].value) ** 2) * + result.params[ + amp_string].stderr)) * 100, 'unit': '%'} amp_string = 's2_amplitude' - result_str_dict['Contrast 2'] = {'value': ((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)*100), - 'error': (np.abs((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)**2 * - result.params['offset'].stderr) + - np.abs((2/(result.params['offset'].value + - result.params[amp_string].value) + (2*result.params[amp_string].value) / - (result.params['offset'].value + result.params[amp_string].value)**2) * - result.params[amp_string].stderr))*100, + result_str_dict['Contrast 2'] = {'value': ((2 * result.params[amp_string].value) / + (result.params['offset'].value + + result.params[amp_string].value) * 100), + 'error': (np.abs( + (2 * result.params[amp_string].value) / + (result.params['offset'].value + result.params[ + amp_string].value) ** 2 * + result.params['offset'].stderr) + + np.abs( + (2 / (result.params['offset'].value + + result.params[ + amp_string].value) + ( + 2 * result.params[ + amp_string].value) / + (result.params['offset'].value + + result.params[ + amp_string].value) ** 2) * + result.params[ + amp_string].stderr)) * 100, 'unit': '%'} - result_str_dict['Phase 1'] = {'value': 180/np.pi*result.params['s1_phase'].value, - 'error': 180/np.pi*result.params['s1_phase'].stderr, - 'unit': 'deg'} + result_str_dict['Phase 1'] = { + 'value': 180 / np.pi * result.params['s1_phase'].value, + 'error': 180 / np.pi * result.params['s1_phase'].stderr, + 'unit': 'deg'} - result_str_dict['Phase 2'] = {'value': 180/np.pi*result.params['s2_phase'].value, - 'error': 180/np.pi*result.params['s2_phase'].stderr, - 'unit': 'deg'} + result_str_dict['Phase 2'] = { + 'value': 180 / np.pi * result.params['s2_phase'].value, + 'error': 180 / np.pi * result.params['s2_phase'].stderr, + 'unit': 'deg'} result_str_dict['Offset'] = {'value': result.params['offset'].value, 'error': result.params['offset'].stderr, @@ -1414,8 +1519,8 @@ def estimate_sinedoublewithexpdecay(self, x_axis, data, params): params['s2_frequency'].set(value=result2.params['frequency'].value) params['s2_phase'].set(value=result2.params['phase'].value) - lifetime = (result1.params['lifetime'].value + result2.params['lifetime'].value)/2 - params['lifetime'].set(value=lifetime, min=2*(x_axis[1]-x_axis[0])) + lifetime = (result1.params['lifetime'].value + result2.params['lifetime'].value) / 2 + params['lifetime'].set(value=lifetime, min=2 * (x_axis[1] - x_axis[0])) params['offset'].set(value=data.mean()) return error, params @@ -1427,7 +1532,8 @@ def estimate_sinedoublewithexpdecay(self, x_axis, data, params): # Problem with stderr: x.stderr will always be 0 for this model! -def make_sinedoublewithtwoexpdecay_fit(self, x_axis, data, estimator, units=None, add_params=None, **kwargs): +def make_sinedoublewithtwoexpdecay_fit(self, x_axis, data, estimator, units=None, + add_params=None, **kwargs): """ Perform a two sine with two exponential decay and offset fit on the provided data. @@ -1451,26 +1557,30 @@ def make_sinedoublewithtwoexpdecay_fit(self, x_axis, data, estimator, units=None params = self._substitute_params(initial_params=params, update_params=add_params) try: - result = two_sine_two_exp_decay_offset.fit(data, x=x_axis, params=params, **kwargs) + result = two_sine_two_exp_decay_offset.fit(data, x=x_axis, params=params, + **kwargs) except: self.log.warning('The sinedoublewithtwoexpdecay fit did not work. ' - 'Error message: {}'.format(str(result.message))) - result = two_sine_two_exp_decay_offset.fit(data, x=x_axis, params=params, **kwargs) + f'Error message: {result.message!s}') + result = two_sine_two_exp_decay_offset.fit(data, x=x_axis, params=params, + **kwargs) if units is None: units = ['arb. unit', 'arb. unit'] - result_str_dict = dict() # create result string for gui or OrderedDict() + result_str_dict = {} # create result string for gui or OrderedDict() period1 = 1 / result.params['e1_frequency'].value try: - period1_err = result.params['e1_frequency'].stderr / (result.params['e1_frequency']) ** 2 + period1_err = result.params['e1_frequency'].stderr / ( + result.params['e1_frequency']) ** 2 except ZeroDivisionError: period1_err = np.inf period2 = 1 / result.params['e2_frequency'].value try: - period2_err = result.params['e2_frequency'].stderr / (result.params['e2_frequency']) ** 2 + period2_err = result.params['e2_frequency'].stderr / ( + result.params['e2_frequency']) ** 2 except ZeroDivisionError: period2_err = np.inf @@ -1484,11 +1594,13 @@ def make_sinedoublewithtwoexpdecay_fit(self, x_axis, data, estimator, units=None result_str_dict['Frequency 1'] = {'value': result.params['e1_frequency'].value, 'error': result.params['e1_frequency'].stderr, - 'unit': 'Hz' if units[0] == 's' else '1/' + units[0]} + 'unit': 'Hz' if units[0] == 's' else '1/' + units[ + 0]} result_str_dict['Frequency 2'] = {'value': result.params['e2_frequency'].value, 'error': result.params['e2_frequency'].stderr, - 'unit': 'Hz' if units[0] == 's' else '1/' + units[0]} + 'unit': 'Hz' if units[0] == 's' else '1/' + units[ + 0]} result_str_dict['Amplitude 1'] = {'value': result.params['e1_amplitude'].value, 'error': result.params['e1_amplitude'].stderr, @@ -1499,36 +1611,58 @@ def make_sinedoublewithtwoexpdecay_fit(self, x_axis, data, estimator, units=None 'unit': units[1]} amp_string = 'e1_amplitude' - result_str_dict['Contrast 1'] = {'value': ((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)*100), - 'error': (np.abs((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)**2 * - result.params['offset'].stderr) + - np.abs((2/(result.params['offset'].value + - result.params[amp_string].value) + (2*result.params[amp_string].value) / - (result.params['offset'].value + result.params[amp_string].value)**2) * - result.params[amp_string].stderr))*100, + result_str_dict['Contrast 1'] = {'value': ((2 * result.params[amp_string].value) / + (result.params['offset'].value + + result.params[amp_string].value) * 100), + 'error': (np.abs( + (2 * result.params[amp_string].value) / + (result.params['offset'].value + result.params[ + amp_string].value) ** 2 * + result.params['offset'].stderr) + + np.abs( + (2 / (result.params['offset'].value + + result.params[ + amp_string].value) + ( + 2 * result.params[ + amp_string].value) / + (result.params['offset'].value + + result.params[ + amp_string].value) ** 2) * + result.params[ + amp_string].stderr)) * 100, 'unit': '%'} amp_string = 'e2_amplitude' - result_str_dict['Contrast 2'] = {'value': ((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)*100), - 'error': (np.abs((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)**2 * - result.params['offset'].stderr) + - np.abs((2/(result.params['offset'].value + - result.params[amp_string].value) + (2*result.params[amp_string].value) / - (result.params['offset'].value + result.params[amp_string].value)**2) * - result.params[amp_string].stderr))*100, + result_str_dict['Contrast 2'] = {'value': ((2 * result.params[amp_string].value) / + (result.params['offset'].value + + result.params[amp_string].value) * 100), + 'error': (np.abs( + (2 * result.params[amp_string].value) / + (result.params['offset'].value + result.params[ + amp_string].value) ** 2 * + result.params['offset'].stderr) + + np.abs( + (2 / (result.params['offset'].value + + result.params[ + amp_string].value) + ( + 2 * result.params[ + amp_string].value) / + (result.params['offset'].value + + result.params[ + amp_string].value) ** 2) * + result.params[ + amp_string].stderr)) * 100, 'unit': '%'} - result_str_dict['Phase 1'] = {'value': 180/np.pi*result.params['e1_phase'].value, - 'error': 180/np.pi*result.params['e1_phase'].stderr, - 'unit': 'deg'} + result_str_dict['Phase 1'] = { + 'value': 180 / np.pi * result.params['e1_phase'].value, + 'error': 180 / np.pi * result.params['e1_phase'].stderr, + 'unit': 'deg'} - result_str_dict['Phase 2'] = {'value': 180/np.pi*result.params['e2_phase'].value, - 'error': 180/np.pi*result.params['e2_phase'].stderr, - 'unit': 'deg'} + result_str_dict['Phase 2'] = { + 'value': 180 / np.pi * result.params['e2_phase'].value, + 'error': 180 / np.pi * result.params['e2_phase'].stderr, + 'unit': 'deg'} result_str_dict['Lifetime 1'] = {'value': result.params['e1_lifetime'].value, 'error': result.params['e1_lifetime'].stderr, @@ -1584,24 +1718,26 @@ def estimate_sinedoublewithtwoexpdecay(self, x_axis, data, params): params['e1_frequency'].set(value=result1.params['frequency'].value) params['e1_phase'].set(value=result1.params['phase'].value) params['e1_lifetime'].set(value=result1.params['lifetime'].value, - min=2*(x_axis[1]-x_axis[0])) + min=2 * (x_axis[1] - x_axis[0])) params['e2_amplitude'].set(value=result2.params['amplitude'].value) params['e2_frequency'].set(value=result2.params['frequency'].value) params['e2_phase'].set(value=result2.params['phase'].value) params['e2_lifetime'].set(value=result2.params['lifetime'].value, - min=2*(x_axis[1]-x_axis[0])) + min=2 * (x_axis[1] - x_axis[0])) params['offset'].set(value=data.mean()) return error, params + ############################################# # Sum of three individual Sinus with offset # ############################################# -def make_sinetriple_fit(self, x_axis, data, estimator, units=None, add_params=None, **kwargs): +def make_sinetriple_fit(self, x_axis, data, estimator, units=None, add_params=None, + **kwargs): """ Perform a three sine with offset fit on the provided data. @param numpy.array x_axis: 1D axis values @@ -1627,29 +1763,32 @@ def make_sinetriple_fit(self, x_axis, data, estimator, units=None, add_params=No result = two_sine_offset.fit(data, x=x_axis, params=params, **kwargs) except: self.log.warning('The threesineexpdecayoffset fit did not work. ' - 'Error message: {}'.format(str(result.message))) + f'Error message: {result.message!s}') result = two_sine_offset.fit(data, x=x_axis, params=params, **kwargs) if units is None: units = ['arb. unit', 'arb. unit'] - result_str_dict = dict() # create result string for gui or OrderedDict() + result_str_dict = {} # create result string for gui or OrderedDict() period1 = 1 / result.params['s1_frequency'].value try: - period1_err = result.params['s1_frequency'].stderr / (result.params['s1_frequency']) ** 2 + period1_err = result.params['s1_frequency'].stderr / ( + result.params['s1_frequency']) ** 2 except ZeroDivisionError: period1_err = np.inf period2 = 1 / result.params['s2_frequency'].value try: - period2_err = result.params['s2_frequency'].stderr / (result.params['s2_frequency']) ** 2 + period2_err = result.params['s2_frequency'].stderr / ( + result.params['s2_frequency']) ** 2 except ZeroDivisionError: period2_err = np.inf period3 = 1 / result.params['s3_frequency'].value try: - period3_err = result.params['s3_frequency'].stderr / (result.params['s3_frequency']) ** 2 + period3_err = result.params['s3_frequency'].stderr / ( + result.params['s3_frequency']) ** 2 except ZeroDivisionError: period3_err = np.inf @@ -1667,15 +1806,18 @@ def make_sinetriple_fit(self, x_axis, data, estimator, units=None, add_params=No result_str_dict['Frequency 1'] = {'value': result.params['s1_frequency'].value, 'error': result.params['s1_frequency'].stderr, - 'unit': 'Hz' if units[0] == 's' else '1/' + units[0]} + 'unit': 'Hz' if units[0] == 's' else '1/' + units[ + 0]} result_str_dict['Frequency 2'] = {'value': result.params['s2_frequency'].value, 'error': result.params['s2_frequency'].stderr, - 'unit': 'Hz' if units[0] == 's' else '1/' + units[0]} + 'unit': 'Hz' if units[0] == 's' else '1/' + units[ + 0]} result_str_dict['Frequency 3'] = {'value': result.params['s3_frequency'].value, 'error': result.params['s3_frequency'].stderr, - 'unit': 'Hz' if units[0] == 's' else '1/' + units[0]} + 'unit': 'Hz' if units[0] == 's' else '1/' + units[ + 0]} result_str_dict['Amplitude 1'] = {'value': result.params['s1_amplitude'].value, 'error': result.params['s1_amplitude'].stderr, @@ -1690,52 +1832,85 @@ def make_sinetriple_fit(self, x_axis, data, estimator, units=None, add_params=No 'unit': units[1]} amp_string = 's1_amplitude' - result_str_dict['Contrast 1'] = {'value': ((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)*100), - 'error': (np.abs((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)**2 * - result.params['offset'].stderr) + - np.abs((2/(result.params['offset'].value + - result.params[amp_string].value) + (2*result.params[amp_string].value) / - (result.params['offset'].value + result.params[amp_string].value)**2) * - result.params[amp_string].stderr))*100, + result_str_dict['Contrast 1'] = {'value': ((2 * result.params[amp_string].value) / + (result.params['offset'].value + + result.params[amp_string].value) * 100), + 'error': (np.abs( + (2 * result.params[amp_string].value) / + (result.params['offset'].value + result.params[ + amp_string].value) ** 2 * + result.params['offset'].stderr) + + np.abs( + (2 / (result.params['offset'].value + + result.params[ + amp_string].value) + ( + 2 * result.params[ + amp_string].value) / + (result.params['offset'].value + + result.params[ + amp_string].value) ** 2) * + result.params[ + amp_string].stderr)) * 100, 'unit': '%'} amp_string = 's2_amplitude' - result_str_dict['Contrast 2'] = {'value': ((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)*100), - 'error': (np.abs((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)**2 * - result.params['offset'].stderr) + - np.abs((2/(result.params['offset'].value + - result.params[amp_string].value) + (2*result.params[amp_string].value) / - (result.params['offset'].value + result.params[amp_string].value)**2) * - result.params[amp_string].stderr))*100, + result_str_dict['Contrast 2'] = {'value': ((2 * result.params[amp_string].value) / + (result.params['offset'].value + + result.params[amp_string].value) * 100), + 'error': (np.abs( + (2 * result.params[amp_string].value) / + (result.params['offset'].value + result.params[ + amp_string].value) ** 2 * + result.params['offset'].stderr) + + np.abs( + (2 / (result.params['offset'].value + + result.params[ + amp_string].value) + ( + 2 * result.params[ + amp_string].value) / + (result.params['offset'].value + + result.params[ + amp_string].value) ** 2) * + result.params[ + amp_string].stderr)) * 100, 'unit': '%'} amp_string = 's3_amplitude' - result_str_dict['Contrast 3'] = {'value': ((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)*100), - 'error': (np.abs((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)**2 * - result.params['offset'].stderr) + - np.abs((2/(result.params['offset'].value + - result.params[amp_string].value) + (2*result.params[amp_string].value) / - (result.params['offset'].value + result.params[amp_string].value)**2) * - result.params[amp_string].stderr))*100, + result_str_dict['Contrast 3'] = {'value': ((2 * result.params[amp_string].value) / + (result.params['offset'].value + + result.params[amp_string].value) * 100), + 'error': (np.abs( + (2 * result.params[amp_string].value) / + (result.params['offset'].value + result.params[ + amp_string].value) ** 2 * + result.params['offset'].stderr) + + np.abs( + (2 / (result.params['offset'].value + + result.params[ + amp_string].value) + ( + 2 * result.params[ + amp_string].value) / + (result.params['offset'].value + + result.params[ + amp_string].value) ** 2) * + result.params[ + amp_string].stderr)) * 100, 'unit': '%'} - result_str_dict['Phase 1'] = {'value': 180/np.pi*result.params['s1_phase'].value, - 'error': 180/np.pi*result.params['s1_phase'].stderr, - 'unit': 'deg'} + result_str_dict['Phase 1'] = { + 'value': 180 / np.pi * result.params['s1_phase'].value, + 'error': 180 / np.pi * result.params['s1_phase'].stderr, + 'unit': 'deg'} - result_str_dict['Phase 2'] = {'value': 180/np.pi*result.params['s2_phase'].value, - 'error': 180/np.pi*result.params['s2_phase'].stderr, - 'unit': 'deg'} + result_str_dict['Phase 2'] = { + 'value': 180 / np.pi * result.params['s2_phase'].value, + 'error': 180 / np.pi * result.params['s2_phase'].stderr, + 'unit': 'deg'} - result_str_dict['Phase 3'] = {'value': 180/np.pi*result.params['s3_phase'].value, - 'error': 180/np.pi*result.params['s3_phase'].stderr, - 'unit': 'deg'} + result_str_dict['Phase 3'] = { + 'value': 180 / np.pi * result.params['s3_phase'].value, + 'error': 180 / np.pi * result.params['s3_phase'].stderr, + 'unit': 'deg'} result_str_dict['Offset'] = {'value': result.params['offset'].value, 'error': result.params['offset'].stderr, @@ -1769,10 +1944,12 @@ def estimate_sinetriple(self, x_axis, data, params): res1 = self.make_sine_fit(x_axis=x_axis, data=data, estimator=self.estimate_sine) data_sub1 = data - res1.best_fit - res2 = self.make_sine_fit(x_axis=x_axis, data=data_sub1, estimator=self.estimate_sine) + res2 = self.make_sine_fit(x_axis=x_axis, data=data_sub1, + estimator=self.estimate_sine) data_sub2 = data_sub1 - res2.best_fit - res3 = self.make_sine_fit(x_axis=x_axis, data=data_sub2, estimator=self.estimate_sine) + res3 = self.make_sine_fit(x_axis=x_axis, data=data_sub2, + estimator=self.estimate_sine) # Fill the parameter dict: params['s1_amplitude'].set(value=res1.params['amplitude'].value) @@ -1791,12 +1968,14 @@ def estimate_sinetriple(self, x_axis, data, params): return error, params + ########################################################################## # Sum of three individual Sinus with offset and single exponential decay # ########################################################################## -def make_sinetriplewithexpdecay_fit(self, x_axis, data, estimator, units=None, add_params=None, **kwargs): +def make_sinetriplewithexpdecay_fit(self, x_axis, data, estimator, units=None, + add_params=None, **kwargs): """ Perform a three sine with one exponential decay offset fit on the provided data. @@ -1819,32 +1998,37 @@ def make_sinetriplewithexpdecay_fit(self, x_axis, data, estimator, units=None, a params = self._substitute_params(initial_params=params, update_params=add_params) try: - result = three_sine_exp_decay_offset.fit(data, x=x_axis, params=params, **kwargs) + result = three_sine_exp_decay_offset.fit(data, x=x_axis, params=params, + **kwargs) except: self.log.warning('The sinetriplewithexpdecay fit did not work. ' - 'Error message: {}'.format(str(result.message))) - result = three_sine_exp_decay_offset.fit(data, x=x_axis, params=params, **kwargs) + f'Error message: {result.message!s}') + result = three_sine_exp_decay_offset.fit(data, x=x_axis, params=params, + **kwargs) if units is None: units = ['arb. unit', 'arb. unit'] - result_str_dict = dict() # create result string for gui or OrderedDict() + result_str_dict = {} # create result string for gui or OrderedDict() period1 = 1 / result.params['s1_frequency'].value try: - period1_err = result.params['s1_frequency'].stderr / (result.params['s1_frequency']) ** 2 + period1_err = result.params['s1_frequency'].stderr / ( + result.params['s1_frequency']) ** 2 except ZeroDivisionError: period1_err = np.inf period2 = 1 / result.params['s2_frequency'].value try: - period2_err = result.params['s2_frequency'].stderr / (result.params['s2_frequency']) ** 2 + period2_err = result.params['s2_frequency'].stderr / ( + result.params['s2_frequency']) ** 2 except ZeroDivisionError: period2_err = np.inf period3 = 1 / result.params['s3_frequency'].value try: - period3_err = result.params['s3_frequency'].stderr / (result.params['s3_frequency']) ** 2 + period3_err = result.params['s3_frequency'].stderr / ( + result.params['s3_frequency']) ** 2 except ZeroDivisionError: period3_err = np.inf @@ -1861,15 +2045,18 @@ def make_sinetriplewithexpdecay_fit(self, x_axis, data, estimator, units=None, a result_str_dict['Frequency 1'] = {'value': result.params['s1_frequency'].value, 'error': result.params['s1_frequency'].stderr, - 'unit': 'Hz' if units[0] == 's' else '1/' + units[0]} + 'unit': 'Hz' if units[0] == 's' else '1/' + units[ + 0]} result_str_dict['Frequency 2'] = {'value': result.params['s2_frequency'].value, 'error': result.params['s2_frequency'].stderr, - 'unit': 'Hz' if units[0] == 's' else '1/' + units[0]} + 'unit': 'Hz' if units[0] == 's' else '1/' + units[ + 0]} result_str_dict['Frequency 3'] = {'value': result.params['s3_frequency'].value, 'error': result.params['s3_frequency'].stderr, - 'unit': 'Hz' if units[0] == 's' else '1/' + units[0]} + 'unit': 'Hz' if units[0] == 's' else '1/' + units[ + 0]} result_str_dict['Amplitude 1'] = {'value': result.params['s1_amplitude'].value, 'error': result.params['s1_amplitude'].stderr, @@ -1883,54 +2070,86 @@ def make_sinetriplewithexpdecay_fit(self, x_axis, data, estimator, units=None, a 'error': result.params['s3_amplitude'].stderr, 'unit': units[1]} - amp_string = 's1_amplitude' - result_str_dict['Contrast 1'] = {'value': ((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)*100), - 'error': (np.abs((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)**2 * - result.params['offset'].stderr) + - np.abs((2/(result.params['offset'].value + - result.params[amp_string].value) + (2*result.params[amp_string].value) / - (result.params['offset'].value + result.params[amp_string].value)**2) * - result.params[amp_string].stderr))*100, + result_str_dict['Contrast 1'] = {'value': ((2 * result.params[amp_string].value) / + (result.params['offset'].value + + result.params[amp_string].value) * 100), + 'error': (np.abs( + (2 * result.params[amp_string].value) / + (result.params['offset'].value + result.params[ + amp_string].value) ** 2 * + result.params['offset'].stderr) + + np.abs( + (2 / (result.params['offset'].value + + result.params[ + amp_string].value) + ( + 2 * result.params[ + amp_string].value) / + (result.params['offset'].value + + result.params[ + amp_string].value) ** 2) * + result.params[ + amp_string].stderr)) * 100, 'unit': '%'} amp_string = 's2_amplitude' - result_str_dict['Contrast 2'] = {'value': ((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)*100), - 'error': (np.abs((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)**2 * - result.params['offset'].stderr) + - np.abs((2/(result.params['offset'].value + - result.params[amp_string].value) + (2*result.params[amp_string].value) / - (result.params['offset'].value + result.params[amp_string].value)**2) * - result.params[amp_string].stderr))*100, + result_str_dict['Contrast 2'] = {'value': ((2 * result.params[amp_string].value) / + (result.params['offset'].value + + result.params[amp_string].value) * 100), + 'error': (np.abs( + (2 * result.params[amp_string].value) / + (result.params['offset'].value + result.params[ + amp_string].value) ** 2 * + result.params['offset'].stderr) + + np.abs( + (2 / (result.params['offset'].value + + result.params[ + amp_string].value) + ( + 2 * result.params[ + amp_string].value) / + (result.params['offset'].value + + result.params[ + amp_string].value) ** 2) * + result.params[ + amp_string].stderr)) * 100, 'unit': '%'} amp_string = 's3_amplitude' - result_str_dict['Contrast 3'] = {'value': ((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)*100), - 'error': (np.abs((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)**2 * - result.params['offset'].stderr) + - np.abs((2/(result.params['offset'].value + - result.params[amp_string].value) + (2*result.params[amp_string].value) / - (result.params['offset'].value + result.params[amp_string].value)**2) * - result.params[amp_string].stderr))*100, + result_str_dict['Contrast 3'] = {'value': ((2 * result.params[amp_string].value) / + (result.params['offset'].value + + result.params[amp_string].value) * 100), + 'error': (np.abs( + (2 * result.params[amp_string].value) / + (result.params['offset'].value + result.params[ + amp_string].value) ** 2 * + result.params['offset'].stderr) + + np.abs( + (2 / (result.params['offset'].value + + result.params[ + amp_string].value) + ( + 2 * result.params[ + amp_string].value) / + (result.params['offset'].value + + result.params[ + amp_string].value) ** 2) * + result.params[ + amp_string].stderr)) * 100, 'unit': '%'} - result_str_dict['Phase 1'] = {'value': 180/np.pi*result.params['s1_phase'].value, - 'error': 180/np.pi*result.params['s1_phase'].stderr, - 'unit': 'deg'} + result_str_dict['Phase 1'] = { + 'value': 180 / np.pi * result.params['s1_phase'].value, + 'error': 180 / np.pi * result.params['s1_phase'].stderr, + 'unit': 'deg'} - result_str_dict['Phase 2'] = {'value': 180/np.pi*result.params['s2_phase'].value, - 'error': 180/np.pi*result.params['s2_phase'].stderr, - 'unit': 'deg'} + result_str_dict['Phase 2'] = { + 'value': 180 / np.pi * result.params['s2_phase'].value, + 'error': 180 / np.pi * result.params['s2_phase'].stderr, + 'unit': 'deg'} - result_str_dict['Phase 3'] = {'value': 180/np.pi*result.params['s3_phase'].value, - 'error': 180/np.pi*result.params['s3_phase'].stderr, - 'unit': 'deg'} + result_str_dict['Phase 3'] = { + 'value': 180 / np.pi * result.params['s3_phase'].value, + 'error': 180 / np.pi * result.params['s3_phase'].stderr, + 'unit': 'deg'} result_str_dict['Lifetime'] = {'value': result.params['lifetime'].value, 'error': result.params['lifetime'].stderr, @@ -1998,19 +2217,22 @@ def estimate_sinetriplewithexpdecay(self, x_axis, data, params): params['s3_frequency'].set(value=res3.params['frequency'].value) params['s3_phase'].set(value=res3.params['phase'].value) - lifetime = (res1.params['lifetime'].value + res2.params['lifetime'].value + res3.params['lifetime'].value)/3 + lifetime = (res1.params['lifetime'].value + res2.params['lifetime'].value + + res3.params['lifetime'].value) / 3 params['lifetime'].set(value=lifetime, - min=2*(x_axis[1]-x_axis[0])) + min=2 * (x_axis[1] - x_axis[0])) params['offset'].set(value=data.mean()) return error, params + ######################################################################### # Sum of three individual Sinus with offset and three exponential decay # ######################################################################### -def make_sinetriplewiththreeexpdecay_fit(self, x_axis, data, estimator, units=None, add_params=None, **kwargs): +def make_sinetriplewiththreeexpdecay_fit(self, x_axis, data, estimator, units=None, + add_params=None, **kwargs): """ Perform a three sine with three exponential decay and offset fit on the provided data. @@ -2034,32 +2256,37 @@ def make_sinetriplewiththreeexpdecay_fit(self, x_axis, data, estimator, units=No params = self._substitute_params(initial_params=params, update_params=add_params) try: - result = three_sine_three_exp_decay_offset.fit(data, x=x_axis, params=params, **kwargs) + result = three_sine_three_exp_decay_offset.fit(data, x=x_axis, params=params, + **kwargs) except: self.log.warning('The twosinetwoexpdecayoffset fit did not work. ' - 'Error message: {}'.format(str(result.message))) - result = three_sine_three_exp_decay_offset.fit(data, x=x_axis, params=params, **kwargs) + f'Error message: {result.message!s}') + result = three_sine_three_exp_decay_offset.fit(data, x=x_axis, params=params, + **kwargs) if units is None: units = ['arb. unit', 'arb. unit'] - result_str_dict = dict() # create result string for gui or OrderedDict() + result_str_dict = {} # create result string for gui or OrderedDict() period1 = 1 / result.params['e1_frequency'].value try: - period1_err = result.params['e1_frequency'].stderr / (result.params['e1_frequency']) ** 2 + period1_err = result.params['e1_frequency'].stderr / ( + result.params['e1_frequency']) ** 2 except ZeroDivisionError: period1_err = np.inf period2 = 1 / result.params['e2_frequency'].value try: - period2_err = result.params['e2_frequency'].stderr / (result.params['e2_frequency']) ** 2 + period2_err = result.params['e2_frequency'].stderr / ( + result.params['e2_frequency']) ** 2 except ZeroDivisionError: period2_err = np.inf period3 = 1 / result.params['e3_frequency'].value try: - period3_err = result.params['e3_frequency'].stderr / (result.params['e3_frequency']) ** 2 + period3_err = result.params['e3_frequency'].stderr / ( + result.params['e3_frequency']) ** 2 except ZeroDivisionError: period3_err = np.inf @@ -2076,15 +2303,18 @@ def make_sinetriplewiththreeexpdecay_fit(self, x_axis, data, estimator, units=No result_str_dict['Frequency 1'] = {'value': result.params['e1_frequency'].value, 'error': result.params['e1_frequency'].stderr, - 'unit': 'Hz' if units[0] == 's' else '1/' + units[0]} + 'unit': 'Hz' if units[0] == 's' else '1/' + units[ + 0]} result_str_dict['Frequency 2'] = {'value': result.params['e2_frequency'].value, 'error': result.params['e2_frequency'].stderr, - 'unit': 'Hz' if units[0] == 's' else '1/' + units[0]} + 'unit': 'Hz' if units[0] == 's' else '1/' + units[ + 0]} result_str_dict['Frequency 3'] = {'value': result.params['e3_frequency'].value, 'error': result.params['e3_frequency'].stderr, - 'unit': 'Hz' if units[0] == 's' else '1/' + units[0]} + 'unit': 'Hz' if units[0] == 's' else '1/' + units[ + 0]} result_str_dict['Amplitude 1'] = {'value': result.params['e1_amplitude'].value, 'error': result.params['e1_amplitude'].stderr, @@ -2099,53 +2329,85 @@ def make_sinetriplewiththreeexpdecay_fit(self, x_axis, data, estimator, units=No 'unit': units[1]} amp_string = 'e1_amplitude' - result_str_dict['Contrast 1'] = {'value': ((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)*100), - 'error': (np.abs((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)**2 * - result.params['offset'].stderr) + - np.abs((2/(result.params['offset'].value + - result.params[amp_string].value) + (2*result.params[amp_string].value) / - (result.params['offset'].value + result.params[amp_string].value)**2) * - result.params[amp_string].stderr))*100, + result_str_dict['Contrast 1'] = {'value': ((2 * result.params[amp_string].value) / + (result.params['offset'].value + + result.params[amp_string].value) * 100), + 'error': (np.abs( + (2 * result.params[amp_string].value) / + (result.params['offset'].value + result.params[ + amp_string].value) ** 2 * + result.params['offset'].stderr) + + np.abs( + (2 / (result.params['offset'].value + + result.params[ + amp_string].value) + ( + 2 * result.params[ + amp_string].value) / + (result.params['offset'].value + + result.params[ + amp_string].value) ** 2) * + result.params[ + amp_string].stderr)) * 100, 'unit': '%'} amp_string = 'e2_amplitude' - result_str_dict['Contrast 2'] = {'value': ((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)*100), - 'error': (np.abs((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)**2 * - result.params['offset'].stderr) + - np.abs((2/(result.params['offset'].value + - result.params[amp_string].value) + (2*result.params[amp_string].value) / - (result.params['offset'].value + result.params[amp_string].value)**2) * - result.params[amp_string].stderr))*100, + result_str_dict['Contrast 2'] = {'value': ((2 * result.params[amp_string].value) / + (result.params['offset'].value + + result.params[amp_string].value) * 100), + 'error': (np.abs( + (2 * result.params[amp_string].value) / + (result.params['offset'].value + result.params[ + amp_string].value) ** 2 * + result.params['offset'].stderr) + + np.abs( + (2 / (result.params['offset'].value + + result.params[ + amp_string].value) + ( + 2 * result.params[ + amp_string].value) / + (result.params['offset'].value + + result.params[ + amp_string].value) ** 2) * + result.params[ + amp_string].stderr)) * 100, 'unit': '%'} amp_string = 'e3_amplitude' - result_str_dict['Contrast 3'] = {'value': ((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)*100), - 'error': (np.abs((2*result.params[amp_string].value) / - (result.params['offset'].value+result.params[amp_string].value)**2 * - result.params['offset'].stderr) + - np.abs((2/(result.params['offset'].value + - result.params[amp_string].value) + (2*result.params[amp_string].value) / - (result.params['offset'].value + result.params[amp_string].value)**2) * - result.params[amp_string].stderr))*100, + result_str_dict['Contrast 3'] = {'value': ((2 * result.params[amp_string].value) / + (result.params['offset'].value + + result.params[amp_string].value) * 100), + 'error': (np.abs( + (2 * result.params[amp_string].value) / + (result.params['offset'].value + result.params[ + amp_string].value) ** 2 * + result.params['offset'].stderr) + + np.abs( + (2 / (result.params['offset'].value + + result.params[ + amp_string].value) + ( + 2 * result.params[ + amp_string].value) / + (result.params['offset'].value + + result.params[ + amp_string].value) ** 2) * + result.params[ + amp_string].stderr)) * 100, 'unit': '%'} + result_str_dict['Phase 1'] = { + 'value': 180 / np.pi * result.params['e1_phase'].value, + 'error': 180 / np.pi * result.params['e1_phase'].stderr, + 'unit': 'deg'} - result_str_dict['Phase 1'] = {'value': 180/np.pi*result.params['e1_phase'].value, - 'error': 180/np.pi*result.params['e1_phase'].stderr, - 'unit': 'deg'} - - result_str_dict['Phase 2'] = {'value': 180/np.pi*result.params['e2_phase'].value, - 'error': 180/np.pi*result.params['e2_phase'].stderr, - 'unit': 'deg'} + result_str_dict['Phase 2'] = { + 'value': 180 / np.pi * result.params['e2_phase'].value, + 'error': 180 / np.pi * result.params['e2_phase'].stderr, + 'unit': 'deg'} - result_str_dict['Phase 3'] = {'value': 180/np.pi*result.params['e3_phase'].value, - 'error': 180/np.pi*result.params['e3_phase'].stderr, - 'unit': 'deg'} + result_str_dict['Phase 3'] = { + 'value': 180 / np.pi * result.params['e3_phase'].value, + 'error': 180 / np.pi * result.params['e3_phase'].stderr, + 'unit': 'deg'} result_str_dict['Lifetime 1'] = {'value': result.params['e1_lifetime'].value, 'error': result.params['e1_lifetime'].stderr, @@ -2211,19 +2473,19 @@ def estimate_sinetriplewiththreeexpdecay(self, x_axis, data, params): params['e1_frequency'].set(value=res1.params['frequency'].value) params['e1_phase'].set(value=res1.params['phase'].value) params['e1_lifetime'].set(value=res1.params['lifetime'].value, - min=2*(x_axis[1]-x_axis[0])) + min=2 * (x_axis[1] - x_axis[0])) params['e2_amplitude'].set(value=res2.params['amplitude'].value) params['e2_frequency'].set(value=res2.params['frequency'].value) params['e2_phase'].set(value=res2.params['phase'].value) params['e2_lifetime'].set(value=res2.params['lifetime'].value, - min=2*(x_axis[1]-x_axis[0])) + min=2 * (x_axis[1] - x_axis[0])) params['e3_amplitude'].set(value=res3.params['amplitude'].value) params['e3_frequency'].set(value=res3.params['frequency'].value) params['e3_phase'].set(value=res3.params['phase'].value) params['e3_lifetime'].set(value=res3.params['lifetime'].value, - min=2*(x_axis[1]-x_axis[0])) + min=2 * (x_axis[1] - x_axis[0])) params['offset'].set(value=data.mean()) diff --git a/qudi_hira_analysis/_qudi_fit_logic.py b/qudi_hira_analysis/_qudi_fit_logic.py index 172b805..7b0acdd 100644 --- a/qudi_hira_analysis/_qudi_fit_logic.py +++ b/qudi_hira_analysis/_qudi_fit_logic.py @@ -30,7 +30,8 @@ import lmfit import numpy as np -logging.basicConfig(format='%(name)s :: %(levelname)s :: %(message)s', level=logging.INFO) +logging.basicConfig(format='%(name)s :: %(levelname)s :: %(message)s', + level=logging.INFO) class FitLogic: @@ -55,26 +56,32 @@ def __init__(self): filenames = [] # for path in directories: - path_list = [os.path.join(os.path.dirname(os.path.realpath(__file__)), '_fitmethods')] + path_list = [ + os.path.join(os.path.dirname(os.path.realpath(__file__)), '_fitmethods')] # adding additional path, to be defined in the config self.log = logging.getLogger(__name__) if self._additional_methods_import_path: if isinstance(self._additional_methods_import_path, str): - self._additional_methods_import_path = [self._additional_methods_import_path] - self.log.info('Adding fit methods path: {}'.format(self._additional_methods_import_path)) + self._additional_methods_import_path = [ + self._additional_methods_import_path] + self.log.info('Adding fit methods path: {}'.format( + self._additional_methods_import_path)) if isinstance(self._additional_methods_import_path, (list, tuple, set)): - self.log.info('Adding fit methods path list: {}'.format(self._additional_methods_import_path)) + self.log.info('Adding fit methods path list: {}'.format( + self._additional_methods_import_path)) for method_import_path in self._additional_methods_import_path: if not os.path.exists(method_import_path): - self.log.error('Specified path "{0}" for import of additional fit methods ' - 'does not exist.'.format(method_import_path)) + self.log.error( + 'Specified path "{}" for import of additional fit methods ' + 'does not exist.'.format(method_import_path)) else: path_list.append(method_import_path) else: - self.log.error('ConfigOption additional_predefined_methods_path needs to either be a string or ' - 'a list of strings.') + self.log.error( + 'ConfigOption additional_predefined_methods_path needs to either be a string or ' + 'a list of strings.') for path in path_list: for f in os.listdir(path): @@ -91,12 +98,12 @@ def __init__(self): # Go through the _fitmethods files and import all methods. # Also determine which methods need to be added to the fit_list dictionary - estimators_for_dict = list() - models_for_dict = list() - fits_for_dict = list() + estimators_for_dict = [] + models_for_dict = [] + fits_for_dict = [] for files in filenames: - mod = importlib.import_module('{0}'.format(files)) + mod = importlib.import_module(f'{files}') for method in dir(mod): ref = getattr(mod, method) if callable(ref) and (inspect.ismethod(ref) or inspect.isfunction(ref)): @@ -105,14 +112,18 @@ def __init__(self): # import methods in Fitlogic setattr(FitLogic, method, ref) # append method to a list of methods to include in the fit_list dictionary - if method_str.startswith('make_') and method_str.endswith('_fit'): - fits_for_dict.append(method_str.split('_', 1)[1].rsplit('_', 1)[0]) - elif method_str.startswith('make_') and method_str.endswith('_model'): - models_for_dict.append(method_str.split('_', 1)[1].rsplit('_', 1)[0]) + if method_str.startswith('make_') and method_str.endswith( + '_fit'): + fits_for_dict.append( + method_str.split('_', 1)[1].rsplit('_', 1)[0]) + elif method_str.startswith('make_') and method_str.endswith( + '_model'): + models_for_dict.append( + method_str.split('_', 1)[1].rsplit('_', 1)[0]) elif method_str.startswith('estimate_'): estimators_for_dict.append(method_str.split('_', 1)[1]) - except: - self.log.error('Method "{0}" could not be imported to FitLogic.' + except Exception as _: + self.log.error('Method "{}" could not be imported to FitLogic.' ''.format(str(method))) fits_for_dict.sort() @@ -138,9 +149,10 @@ def __init__(self): # Attach make_*_model method to fit_list if fit_name in models_for_dict: - self.fit_list[dimension][fit_name]['make_model'] = getattr(self, model_method) + self.fit_list[dimension][fit_name]['make_model'] = getattr(self, + model_method) else: - self.log.error('No make_*_model method for fit "{0}" found in FitLogic.' + self.log.error('No make_*_model method for fit "{}" found in FitLogic.' ''.format(fit_name)) # Attach all estimate_* methods to corresponding fit method in fit_list @@ -148,14 +160,16 @@ def __init__(self): for estimator_name in estimators_for_dict: estimator_method = 'estimate_' + estimator_name if fit_name == estimator_name: - self.fit_list[dimension][fit_name]['generic'] = getattr(self, estimator_method) + self.fit_list[dimension][fit_name]['generic'] = getattr(self, + estimator_method) found_estimator = True elif estimator_name.startswith(fit_name + '_'): custom_name = estimator_name.split('_', 1)[1] - self.fit_list[dimension][fit_name][custom_name] = getattr(self, estimator_method) + self.fit_list[dimension][fit_name][custom_name] = getattr(self, + estimator_method) found_estimator = True if not found_estimator: - self.log.error('No estimator method for fit "{0}" found in FitLogic.' + self.log.error('No estimator method for fit "{}" found in FitLogic.' ''.format(fit_name)) # self.log.info('Methods were included to FitLogic, but only if naming is right: check the' @@ -171,7 +185,6 @@ def on_activate(self): def on_deactivate(self): """ """ - pass def validate_load_fits(self, fits): """ Take fit names and estimators from a dict and check if they are valid. @@ -217,7 +230,7 @@ def validate_load_fits(self, fits): new_fit['parameters'] = par user_fits[dim][name] = new_fit except KeyError: - self.log.exception('Failed to validate fit {0}'.format(name)) + self.log.exception(f'Failed to validate fit {name}') continue return user_fits @@ -236,11 +249,13 @@ def prepare_save_fits(self, fits): save_fits[dim] = OrderedDict() for name, fit in dfits.items(): try: - new_fit = {'fit_function': fit['fit_name'], 'estimator': fit['est_name'], + new_fit = {'fit_function': fit['fit_name'], + 'estimator': fit['est_name'], 'parameters': fit['parameters'].dumps()} save_fits[dim][name] = new_fit except KeyError: - self.log.exception('Error while preparing fit {0} for saving.'.format(name)) + self.log.exception( + f'Error while preparing fit {name} for saving.') continue return save_fits @@ -258,7 +273,7 @@ def make_fit_container(self, container_name, dimension): return FitContainer(self, container_name, dimension) -class FitContainer(): +class FitContainer: """ A class for managing a single flexible fit setting in a logic module. """ @@ -280,7 +295,7 @@ def __init__(self, fit_logic, name, dimension): elif dimension == '3d': self.dim = 3 else: - raise Exception('Invalid dimension {0}'.format(dimension)) + raise Exception(f'Invalid dimension {dimension}') self.dimension = dimension self.fit_list = OrderedDict() # variables for fitting @@ -289,7 +304,7 @@ def __init__(self, fit_logic, name, dimension): self.current_fit_param = lmfit.parameter.Parameters() self.current_fit_result = None self.use_settings = None - self.units = ['independent variable {0}'.format(i + 1) for i in range(self.dim)] + self.units = [f'independent variable {i + 1}' for i in range(self.dim)] self.units.append('dependent variable') def set_units(self, units): @@ -340,7 +355,8 @@ def set_current_fit(self, current_fit): This is a reserved name that will do nothing and should not display a fit line if set. """ if current_fit not in self.fit_list and current_fit != 'No Fit': - self.fit_logic.log.warning('{0} not in {1} fit list!'.format(current_fit, self.name)) + self.fit_logic.log.warning( + f'{current_fit} not in {self.name} fit list!') self.current_fit = 'No Fit' else: self.current_fit = current_fit @@ -350,13 +366,18 @@ def set_current_fit(self, current_fit): # Update the use parameter dictionary for para in use_settings: if use_settings[para]: - self.use_settings[para] = self.fit_list[self.current_fit]['parameters'][para] + self.use_settings[para] = \ + self.fit_list[self.current_fit]['parameters'][para] else: self.use_settings = None self.clear_result() return self.current_fit, self.use_settings - def do_fit(self, x_data, y_data): + def do_fit( + self, + x_data: np.ndarray, + y_data: np.ndarray + ) -> tuple[np.ndarray, np.ndarray, lmfit.model.ModelResult]: """Performs the chosen fit on the measured data. @param array x_data: optional, 1D np.array or 1D list with the x values. If None is passed then the module x values are @@ -388,6 +409,7 @@ def do_fit(self, x_data, y_data): start=x_data[0], stop=x_data[-1], num=int(len(x_data) * self.fit_granularity_fact)) + fit_y = None # set the keyword arguments, which will be passed to the fit. kwargs = { @@ -408,7 +430,7 @@ def do_fit(self, x_data, y_data): else: self.fit_logic.log.warning( - 'The Fit Function "{0}" is not implemented to be used in the ODMR Logic. ' + 'The Fit Function "{}" is not implemented to be used in the ODMR Logic. ' 'Correct that! Fit Call will be skipped and Fit Function will be set to ' '"No Fit".'.format(self.current_fit))