summaryrefslogtreecommitdiff
path: root/silx/math/fit/fittheories.py
diff options
context:
space:
mode:
Diffstat (limited to 'silx/math/fit/fittheories.py')
-rw-r--r--silx/math/fit/fittheories.py1374
1 files changed, 0 insertions, 1374 deletions
diff --git a/silx/math/fit/fittheories.py b/silx/math/fit/fittheories.py
deleted file mode 100644
index 6b19a38..0000000
--- a/silx/math/fit/fittheories.py
+++ /dev/null
@@ -1,1374 +0,0 @@
-# coding: utf-8
-#/*##########################################################################
-#
-# Copyright (c) 2004-2020 European Synchrotron Radiation Facility
-#
-# Permission is hereby granted, free of charge, to any person obtaining a copy
-# of this software and associated documentation files (the "Software"), to deal
-# in the Software without restriction, including without limitation the rights
-# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
-# copies of the Software, and to permit persons to whom the Software is
-# furnished to do so, subject to the following conditions:
-#
-# The above copyright notice and this permission notice shall be included in
-# all copies or substantial portions of the Software.
-#
-# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
-# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
-# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
-# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
-# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
-# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
-# THE SOFTWARE.
-#
-########################################################################### */
-"""This modules provides a set of fit functions and associated
-estimation functions in a format that can be imported into a
-:class:`silx.math.fit.FitManager` instance.
-
-These functions are well suited for fitting multiple gaussian shaped peaks
-typically found in spectroscopy data. The estimation functions are designed
-to detect how many peaks are present in the data, and provide an initial
-estimate for their height, their center location and their full-width
-at half maximum (fwhm).
-
-The limitation of these estimation algorithms is that only gaussians having a
-similar fwhm can be detected by the peak search algorithm.
-This *search fwhm* can be defined by the user, if
-he knows the characteristics of his data, or can be automatically estimated
-based on the fwhm of the largest peak in the data.
-
-The source code of this module can serve as template for defining your own
-fit functions.
-
-The functions to be imported by :meth:`FitManager.loadtheories` are defined by
-a dictionary :const:`THEORY`: with the following structure::
-
- from silx.math.fit.fittheory import FitTheory
-
- THEORY = {
- 'theory_name_1': FitTheory(
- description='Description of theory 1',
- function=fitfunction1,
- parameters=('param name 1', 'param name 2', …),
- estimate=estimation_function1,
- configure=configuration_function1,
- derivative=derivative_function1),
-
- 'theory_name_2': FitTheory(…),
- }
-
-.. note::
-
- Consider using an OrderedDict instead of a regular dictionary, when
- defining your own theory dictionary, if the order matters to you.
- This will likely be the case if you intend to load a selection of
- functions in a GUI such as :class:`silx.gui.fit.FitManager`.
-
-Theory names can be customized (e.g. ``gauss, lorentz, splitgauss``…).
-
-The mandatory parameters for :class:`FitTheory` are ``function`` and
-``parameters``.
-
-You can also define an ``INIT`` function that will be executed by
-:meth:`FitManager.loadtheories`.
-
-See the documentation of :class:`silx.math.fit.fittheory.FitTheory`
-for more information.
-
-Module members:
----------------
-"""
-import numpy
-from collections import OrderedDict
-import logging
-
-from silx.math.fit import functions
-from silx.math.fit.peaks import peak_search, guess_fwhm
-from silx.math.fit.filters import strip, savitsky_golay
-from silx.math.fit.leastsq import leastsq
-from silx.math.fit.fittheory import FitTheory
-
-_logger = logging.getLogger(__name__)
-
-__authors__ = ["V.A. Sole", "P. Knobel"]
-__license__ = "MIT"
-__date__ = "15/05/2017"
-
-
-DEFAULT_CONFIG = {
- 'NoConstraintsFlag': False,
- 'PositiveFwhmFlag': True,
- 'PositiveHeightAreaFlag': True,
- 'SameFwhmFlag': False,
- 'QuotedPositionFlag': False, # peak not outside data range
- 'QuotedEtaFlag': False, # force 0 < eta < 1
- # Peak detection
- 'AutoScaling': False,
- 'Yscaling': 1.0,
- 'FwhmPoints': 8,
- 'AutoFwhm': True,
- 'Sensitivity': 2.5,
- 'ForcePeakPresence': True,
- # Hypermet
- 'HypermetTails': 15,
- 'QuotedFwhmFlag': 0,
- 'MaxFwhm2InputRatio': 1.5,
- 'MinFwhm2InputRatio': 0.4,
- # short tail parameters
- 'MinGaussArea4ShortTail': 50000.,
- 'InitialShortTailAreaRatio': 0.050,
- 'MaxShortTailAreaRatio': 0.100,
- 'MinShortTailAreaRatio': 0.0010,
- 'InitialShortTailSlopeRatio': 0.70,
- 'MaxShortTailSlopeRatio': 2.00,
- 'MinShortTailSlopeRatio': 0.50,
- # long tail parameters
- 'MinGaussArea4LongTail': 1000.0,
- 'InitialLongTailAreaRatio': 0.050,
- 'MaxLongTailAreaRatio': 0.300,
- 'MinLongTailAreaRatio': 0.010,
- 'InitialLongTailSlopeRatio': 20.0,
- 'MaxLongTailSlopeRatio': 50.0,
- 'MinLongTailSlopeRatio': 5.0,
- # step tail
- 'MinGaussHeight4StepTail': 5000.,
- 'InitialStepTailHeightRatio': 0.002,
- 'MaxStepTailHeightRatio': 0.0100,
- 'MinStepTailHeightRatio': 0.0001,
- # Hypermet constraints
- # position in range [estimated position +- estimated fwhm/2]
- 'HypermetQuotedPositionFlag': True,
- 'DeltaPositionFwhmUnits': 0.5,
- 'SameSlopeRatioFlag': 1,
- 'SameAreaRatioFlag': 1,
- # Strip bg removal
- 'StripBackgroundFlag': True,
- 'SmoothingFlag': True,
- 'SmoothingWidth': 5,
- 'StripWidth': 2,
- 'StripIterations': 5000,
- 'StripThresholdFactor': 1.0}
-"""This dictionary defines default configuration parameters that have effects
-on fit functions and estimation functions, mainly on fit constraints.
-This dictionary is accessible as attribute :attr:`FitTheories.config`,
-which can be modified by configuration functions defined in
-:const:`CONFIGURE`.
-"""
-
-CFREE = 0
-CPOSITIVE = 1
-CQUOTED = 2
-CFIXED = 3
-CFACTOR = 4
-CDELTA = 5
-CSUM = 6
-CIGNORED = 7
-
-
-class FitTheories(object):
- """Class wrapping functions from :class:`silx.math.fit.functions`
- and providing estimate functions for all of these fit functions."""
- def __init__(self, config=None):
- if config is None:
- self.config = DEFAULT_CONFIG
- else:
- self.config = config
-
- def ahypermet(self, x, *pars):
- """
- Wrapping of :func:`silx.math.fit.functions.sum_ahypermet` without
- the tail flags in the function signature.
-
- Depending on the value of `self.config['HypermetTails']`, one can
- activate or deactivate the various terms of the hypermet function.
-
- `self.config['HypermetTails']` must be an integer between 0 and 15.
- It is a set of 4 binary flags, one for activating each one of the
- hypermet terms: *gaussian function, short tail, long tail, step*.
-
- For example, 15 can be expressed as ``1111`` in base 2, so a flag of
- 15 means all terms are active.
- """
- g_term = self.config['HypermetTails'] & 1
- st_term = (self.config['HypermetTails'] >> 1) & 1
- lt_term = (self.config['HypermetTails'] >> 2) & 1
- step_term = (self.config['HypermetTails'] >> 3) & 1
- return functions.sum_ahypermet(x, *pars,
- gaussian_term=g_term, st_term=st_term,
- lt_term=lt_term, step_term=step_term)
-
- def poly(self, x, *pars):
- """Order n polynomial.
- The order of the polynomial is defined by the number of
- coefficients (``*pars``).
-
- """
- p = numpy.poly1d(pars)
- return p(x)
-
- @staticmethod
- def estimate_poly(x, y, n=2):
- """Estimate polynomial coefficients for a degree n polynomial.
-
- """
- pcoeffs = numpy.polyfit(x, y, n)
- constraints = numpy.zeros((n + 1, 3), numpy.float64)
- return pcoeffs, constraints
-
- def estimate_quadratic(self, x, y):
- """Estimate quadratic coefficients
-
- """
- return self.estimate_poly(x, y, n=2)
-
- def estimate_cubic(self, x, y):
- """Estimate coefficients for a degree 3 polynomial
-
- """
- return self.estimate_poly(x, y, n=3)
-
- def estimate_quartic(self, x, y):
- """Estimate coefficients for a degree 4 polynomial
-
- """
- return self.estimate_poly(x, y, n=4)
-
- def estimate_quintic(self, x, y):
- """Estimate coefficients for a degree 5 polynomial
-
- """
- return self.estimate_poly(x, y, n=5)
-
- def strip_bg(self, y):
- """Return the strip background of y, using parameters from
- :attr:`config` dictionary (*StripBackgroundFlag, StripWidth,
- StripIterations, StripThresholdFactor*)"""
- remove_strip_bg = self.config.get('StripBackgroundFlag', False)
- if remove_strip_bg:
- if self.config['SmoothingFlag']:
- y = savitsky_golay(y, self.config['SmoothingWidth'])
- strip_width = self.config['StripWidth']
- strip_niterations = self.config['StripIterations']
- strip_thr_factor = self.config['StripThresholdFactor']
- return strip(y, w=strip_width,
- niterations=strip_niterations,
- factor=strip_thr_factor)
- else:
- return numpy.zeros_like(y)
-
- def guess_yscaling(self, y):
- """Estimate scaling for y prior to peak search.
- A smoothing filter is applied to y to estimate the noise level
- (chi-squared)
-
- :param y: Data array
- :return: Scaling factor
- """
- # ensure y is an array
- yy = numpy.array(y, copy=False)
-
- # smooth
- convolution_kernel = numpy.ones(shape=(3,)) / 3.
- ysmooth = numpy.convolve(y, convolution_kernel, mode="same")
-
- # remove zeros
- idx_array = numpy.fabs(y) > 0.0
- yy = yy[idx_array]
- ysmooth = ysmooth[idx_array]
-
- # compute scaling factor
- chisq = numpy.mean((yy - ysmooth)**2 / numpy.fabs(yy))
- if chisq > 0:
- return 1. / chisq
- else:
- return 1.0
-
- def peak_search(self, y, fwhm, sensitivity):
- """Search for peaks in y array, after padding the array and
- multiplying its value by a scaling factor.
-
- :param y: 1-D data array
- :param int fwhm: Typical full width at half maximum for peaks,
- in number of points. This parameter is used for to discriminate between
- true peaks and background fluctuations.
- :param float sensitivity: Sensitivity parameter. This is a threshold factor
- for peak detection. Only peaks larger than the standard deviation
- of the noise multiplied by this sensitivity parameter are detected.
- :return: List of peak indices
- """
- # add padding
- ysearch = numpy.ones((len(y) + 2 * fwhm,), numpy.float64)
- ysearch[0:fwhm] = y[0]
- ysearch[-1:-fwhm - 1:-1] = y[len(y)-1]
- ysearch[fwhm:fwhm + len(y)] = y[:]
-
- scaling = self.guess_yscaling(y) if self.config["AutoScaling"] else self.config["Yscaling"]
-
- if len(ysearch) > 1.5 * fwhm:
- peaks = peak_search(scaling * ysearch,
- fwhm=fwhm, sensitivity=sensitivity)
- return [peak_index - fwhm for peak_index in peaks
- if 0 <= peak_index - fwhm < len(y)]
- else:
- return []
-
- def estimate_height_position_fwhm(self, x, y):
- """Estimation of *Height, Position, FWHM* of peaks, for gaussian-like
- curves.
-
- This functions finds how many parameters are needed, based on the
- number of peaks detected. Then it estimates the fit parameters
- with a few iterations of fitting gaussian functions.
-
- :param x: Array of abscissa values
- :param y: Array of ordinate values (``y = f(x)``)
- :return: Tuple of estimated fit parameters and fit constraints.
- Parameters to be estimated for each peak are:
- *Height, Position, FWHM*.
- Fit constraints depend on :attr:`config`.
- """
- fittedpar = []
-
- bg = self.strip_bg(y)
-
- if self.config['AutoFwhm']:
- search_fwhm = guess_fwhm(y)
- else:
- search_fwhm = int(float(self.config['FwhmPoints']))
- search_sens = float(self.config['Sensitivity'])
-
- if search_fwhm < 3:
- _logger.warning("Setting peak fwhm to 3 (lower limit)")
- search_fwhm = 3
- self.config['FwhmPoints'] = 3
-
- if search_sens < 1:
- _logger.warning("Setting peak search sensitivity to 1. " +
- "(lower limit to filter out noise peaks)")
- search_sens = 1
- self.config['Sensitivity'] = 1
-
- npoints = len(y)
-
- # Find indices of peaks in data array
- peaks = self.peak_search(y,
- fwhm=search_fwhm,
- sensitivity=search_sens)
-
- if not len(peaks):
- forcepeak = int(float(self.config.get('ForcePeakPresence', 0)))
- if forcepeak:
- delta = y - bg
- # get index of global maximum
- # (first one if several samples are equal to this value)
- peaks = [numpy.nonzero(delta == delta.max())[0][0]]
-
- # Find index of largest peak in peaks array
- index_largest_peak = 0
- if len(peaks) > 0:
- # estimate fwhm as 5 * sampling interval
- sig = 5 * abs(x[npoints - 1] - x[0]) / npoints
- peakpos = x[int(peaks[0])]
- if abs(peakpos) < 1.0e-16:
- peakpos = 0.0
- param = numpy.array(
- [y[int(peaks[0])] - bg[int(peaks[0])], peakpos, sig])
- height_largest_peak = param[0]
- peak_index = 1
- for i in peaks[1:]:
- param2 = numpy.array(
- [y[int(i)] - bg[int(i)], x[int(i)], sig])
- param = numpy.concatenate((param, param2))
- if param2[0] > height_largest_peak:
- height_largest_peak = param2[0]
- index_largest_peak = peak_index
- peak_index += 1
-
- # Subtract background
- xw = x
- yw = y - bg
-
- cons = numpy.zeros((len(param), 3), numpy.float64)
-
- # peak height must be positive
- cons[0:len(param):3, 0] = CPOSITIVE
- # force peaks to stay around their position
- cons[1:len(param):3, 0] = CQUOTED
-
- # set possible peak range to estimated peak +- guessed fwhm
- if len(xw) > search_fwhm:
- fwhmx = numpy.fabs(xw[int(search_fwhm)] - xw[0])
- cons[1:len(param):3, 1] = param[1:len(param):3] - 0.5 * fwhmx
- cons[1:len(param):3, 2] = param[1:len(param):3] + 0.5 * fwhmx
- else:
- shape = [max(1, int(x)) for x in (param[1:len(param):3])]
- cons[1:len(param):3, 1] = min(xw) * numpy.ones(
- shape,
- numpy.float64)
- cons[1:len(param):3, 2] = max(xw) * numpy.ones(
- shape,
- numpy.float64)
-
- # ensure fwhm is positive
- cons[2:len(param):3, 0] = CPOSITIVE
-
- # run a quick iterative fit (4 iterations) to improve
- # estimations
- fittedpar, _, _ = leastsq(functions.sum_gauss, xw, yw, param,
- max_iter=4, constraints=cons.tolist(),
- full_output=True)
-
- # set final constraints based on config parameters
- cons = numpy.zeros((len(fittedpar), 3), numpy.float64)
- peak_index = 0
- for i in range(len(peaks)):
- # Setup height area constrains
- if not self.config['NoConstraintsFlag']:
- if self.config['PositiveHeightAreaFlag']:
- cons[peak_index, 0] = CPOSITIVE
- cons[peak_index, 1] = 0
- cons[peak_index, 2] = 0
- peak_index += 1
-
- # Setup position constrains
- if not self.config['NoConstraintsFlag']:
- if self.config['QuotedPositionFlag']:
- cons[peak_index, 0] = CQUOTED
- cons[peak_index, 1] = min(x)
- cons[peak_index, 2] = max(x)
- peak_index += 1
-
- # Setup positive FWHM constrains
- if not self.config['NoConstraintsFlag']:
- if self.config['PositiveFwhmFlag']:
- cons[peak_index, 0] = CPOSITIVE
- cons[peak_index, 1] = 0
- cons[peak_index, 2] = 0
- if self.config['SameFwhmFlag']:
- if i != index_largest_peak:
- cons[peak_index, 0] = CFACTOR
- cons[peak_index, 1] = 3 * index_largest_peak + 2
- cons[peak_index, 2] = 1.0
- peak_index += 1
-
- return fittedpar, cons
-
- def estimate_agauss(self, x, y):
- """Estimation of *Area, Position, FWHM* of peaks, for gaussian-like
- curves.
-
- This functions uses :meth:`estimate_height_position_fwhm`, then
- converts the height parameters to area under the curve with the
- formula ``area = sqrt(2*pi) * height * fwhm / (2 * sqrt(2 * log(2))``
-
- :param x: Array of abscissa values
- :param y: Array of ordinate values (``y = f(x)``)
- :return: Tuple of estimated fit parameters and fit constraints.
- Parameters to be estimated for each peak are:
- *Area, Position, FWHM*.
- Fit constraints depend on :attr:`config`.
- """
- fittedpar, cons = self.estimate_height_position_fwhm(x, y)
- # get the number of found peaks
- npeaks = len(fittedpar) // 3
- for i in range(npeaks):
- height = fittedpar[3 * i]
- fwhm = fittedpar[3 * i + 2]
- # Replace height with area in fittedpar
- fittedpar[3 * i] = numpy.sqrt(2 * numpy.pi) * height * fwhm / (
- 2.0 * numpy.sqrt(2 * numpy.log(2)))
- return fittedpar, cons
-
- def estimate_alorentz(self, x, y):
- """Estimation of *Area, Position, FWHM* of peaks, for Lorentzian
- curves.
-
- This functions uses :meth:`estimate_height_position_fwhm`, then
- converts the height parameters to area under the curve with the
- formula ``area = height * fwhm * 0.5 * pi``
-
- :param x: Array of abscissa values
- :param y: Array of ordinate values (``y = f(x)``)
- :return: Tuple of estimated fit parameters and fit constraints.
- Parameters to be estimated for each peak are:
- *Area, Position, FWHM*.
- Fit constraints depend on :attr:`config`.
- """
- fittedpar, cons = self.estimate_height_position_fwhm(x, y)
- # get the number of found peaks
- npeaks = len(fittedpar) // 3
- for i in range(npeaks):
- height = fittedpar[3 * i]
- fwhm = fittedpar[3 * i + 2]
- # Replace height with area in fittedpar
- fittedpar[3 * i] = (height * fwhm * 0.5 * numpy.pi)
- return fittedpar, cons
-
- def estimate_splitgauss(self, x, y):
- """Estimation of *Height, Position, FWHM1, FWHM2* of peaks, for
- asymmetric gaussian-like curves.
-
- This functions uses :meth:`estimate_height_position_fwhm`, then
- adds a second (identical) estimation of FWHM to the fit parameters
- for each peak, and the corresponding constraint.
-
- :param x: Array of abscissa values
- :param y: Array of ordinate values (``y = f(x)``)
- :return: Tuple of estimated fit parameters and fit constraints.
- Parameters to be estimated for each peak are:
- *Height, Position, FWHM1, FWHM2*.
- Fit constraints depend on :attr:`config`.
- """
- fittedpar, cons = self.estimate_height_position_fwhm(x, y)
- # get the number of found peaks
- npeaks = len(fittedpar) // 3
- estimated_parameters = []
- estimated_constraints = numpy.zeros((4 * npeaks, 3), numpy.float64)
- for i in range(npeaks):
- for j in range(3):
- estimated_parameters.append(fittedpar[3 * i + j])
- # fwhm2 estimate = fwhm1
- estimated_parameters.append(fittedpar[3 * i + 2])
- # height
- estimated_constraints[4 * i, 0] = cons[3 * i, 0]
- estimated_constraints[4 * i, 1] = cons[3 * i, 1]
- estimated_constraints[4 * i, 2] = cons[3 * i, 2]
- # position
- estimated_constraints[4 * i + 1, 0] = cons[3 * i + 1, 0]
- estimated_constraints[4 * i + 1, 1] = cons[3 * i + 1, 1]
- estimated_constraints[4 * i + 1, 2] = cons[3 * i + 1, 2]
- # fwhm1
- estimated_constraints[4 * i + 2, 0] = cons[3 * i + 2, 0]
- estimated_constraints[4 * i + 2, 1] = cons[3 * i + 2, 1]
- estimated_constraints[4 * i + 2, 2] = cons[3 * i + 2, 2]
- # fwhm2
- estimated_constraints[4 * i + 3, 0] = cons[3 * i + 2, 0]
- estimated_constraints[4 * i + 3, 1] = cons[3 * i + 2, 1]
- estimated_constraints[4 * i + 3, 2] = cons[3 * i + 2, 2]
- if cons[3 * i + 2, 0] == CFACTOR:
- # convert indices of related parameters
- # (this happens if SameFwhmFlag == True)
- estimated_constraints[4 * i + 2, 1] = \
- int(cons[3 * i + 2, 1] / 3) * 4 + 2
- estimated_constraints[4 * i + 3, 1] = \
- int(cons[3 * i + 2, 1] / 3) * 4 + 3
- return estimated_parameters, estimated_constraints
-
- def estimate_pvoigt(self, x, y):
- """Estimation of *Height, Position, FWHM, eta* of peaks, for
- pseudo-Voigt curves.
-
- Pseudo-Voigt are a sum of a gaussian curve *G(x)* and a lorentzian
- curve *L(x)* with the same height, center, fwhm parameters:
- ``y(x) = eta * G(x) + (1-eta) * L(x)``
-
- This functions uses :meth:`estimate_height_position_fwhm`, then
- adds a constant estimation of *eta* (0.5) to the fit parameters
- for each peak, and the corresponding constraint.
-
- :param x: Array of abscissa values
- :param y: Array of ordinate values (``y = f(x)``)
- :return: Tuple of estimated fit parameters and fit constraints.
- Parameters to be estimated for each peak are:
- *Height, Position, FWHM, eta*.
- Constraint for the eta parameter can be set to QUOTED (0.--1.)
- by setting :attr:`config`['QuotedEtaFlag'] to ``True``.
- If this is not the case, the constraint code is set to FREE.
- """
- fittedpar, cons = self.estimate_height_position_fwhm(x, y)
- npeaks = len(fittedpar) // 3
- newpar = []
- newcons = numpy.zeros((4 * npeaks, 3), numpy.float64)
- # find out related parameters proper index
- if not self.config['NoConstraintsFlag']:
- if self.config['SameFwhmFlag']:
- j = 0
- # get the index of the free FWHM
- for i in range(npeaks):
- if cons[3 * i + 2, 0] != 4:
- j = i
- for i in range(npeaks):
- if i != j:
- cons[3 * i + 2, 1] = 4 * j + 2
- for i in range(npeaks):
- newpar.append(fittedpar[3 * i])
- newpar.append(fittedpar[3 * i + 1])
- newpar.append(fittedpar[3 * i + 2])
- newpar.append(0.5)
- # height
- newcons[4 * i, 0] = cons[3 * i, 0]
- newcons[4 * i, 1] = cons[3 * i, 1]
- newcons[4 * i, 2] = cons[3 * i, 2]
- # position
- newcons[4 * i + 1, 0] = cons[3 * i + 1, 0]
- newcons[4 * i + 1, 1] = cons[3 * i + 1, 1]
- newcons[4 * i + 1, 2] = cons[3 * i + 1, 2]
- # fwhm
- newcons[4 * i + 2, 0] = cons[3 * i + 2, 0]
- newcons[4 * i + 2, 1] = cons[3 * i + 2, 1]
- newcons[4 * i + 2, 2] = cons[3 * i + 2, 2]
- # Eta constrains
- newcons[4 * i + 3, 0] = CFREE
- newcons[4 * i + 3, 1] = 0
- newcons[4 * i + 3, 2] = 0
- if self.config['QuotedEtaFlag']:
- newcons[4 * i + 3, 0] = CQUOTED
- newcons[4 * i + 3, 1] = 0.0
- newcons[4 * i + 3, 2] = 1.0
- return newpar, newcons
-
- def estimate_splitpvoigt(self, x, y):
- """Estimation of *Height, Position, FWHM1, FWHM2, eta* of peaks, for
- asymmetric pseudo-Voigt curves.
-
- This functions uses :meth:`estimate_height_position_fwhm`, then
- adds an identical FWHM2 parameter and a constant estimation of
- *eta* (0.5) to the fit parameters for each peak, and the corresponding
- constraints.
-
- Constraint for the eta parameter can be set to QUOTED (0.--1.)
- by setting :attr:`config`['QuotedEtaFlag'] to ``True``.
- If this is not the case, the constraint code is set to FREE.
-
- :param x: Array of abscissa values
- :param y: Array of ordinate values (``y = f(x)``)
- :return: Tuple of estimated fit parameters and fit constraints.
- Parameters to be estimated for each peak are:
- *Height, Position, FWHM1, FWHM2, eta*.
- """
- fittedpar, cons = self.estimate_height_position_fwhm(x, y)
- npeaks = len(fittedpar) // 3
- newpar = []
- newcons = numpy.zeros((5 * npeaks, 3), numpy.float64)
- # find out related parameters proper index
- if not self.config['NoConstraintsFlag']:
- if self.config['SameFwhmFlag']:
- j = 0
- # get the index of the free FWHM
- for i in range(npeaks):
- if cons[3 * i + 2, 0] != 4:
- j = i
- for i in range(npeaks):
- if i != j:
- cons[3 * i + 2, 1] = 4 * j + 2
- for i in range(npeaks):
- # height
- newpar.append(fittedpar[3 * i])
- # position
- newpar.append(fittedpar[3 * i + 1])
- # fwhm1
- newpar.append(fittedpar[3 * i + 2])
- # fwhm2 estimate equal to fwhm1
- newpar.append(fittedpar[3 * i + 2])
- # eta
- newpar.append(0.5)
- # constraint codes
- # ----------------
- # height
- newcons[5 * i, 0] = cons[3 * i, 0]
- # position
- newcons[5 * i + 1, 0] = cons[3 * i + 1, 0]
- # fwhm1
- newcons[5 * i + 2, 0] = cons[3 * i + 2, 0]
- # fwhm2
- newcons[5 * i + 3, 0] = cons[3 * i + 2, 0]
- # cons 1
- # ------
- newcons[5 * i, 1] = cons[3 * i, 1]
- newcons[5 * i + 1, 1] = cons[3 * i + 1, 1]
- newcons[5 * i + 2, 1] = cons[3 * i + 2, 1]
- newcons[5 * i + 3, 1] = cons[3 * i + 2, 1]
- # cons 2
- # ------
- newcons[5 * i, 2] = cons[3 * i, 2]
- newcons[5 * i + 1, 2] = cons[3 * i + 1, 2]
- newcons[5 * i + 2, 2] = cons[3 * i + 2, 2]
- newcons[5 * i + 3, 2] = cons[3 * i + 2, 2]
-
- if cons[3 * i + 2, 0] == CFACTOR:
- # fwhm2 connstraint depends on fwhm1
- newcons[5 * i + 3, 1] = newcons[5 * i + 2, 1] + 1
- # eta constraints
- newcons[5 * i + 4, 0] = CFREE
- newcons[5 * i + 4, 1] = 0
- newcons[5 * i + 4, 2] = 0
- if self.config['QuotedEtaFlag']:
- newcons[5 * i + 4, 0] = CQUOTED
- newcons[5 * i + 4, 1] = 0.0
- newcons[5 * i + 4, 2] = 1.0
- return newpar, newcons
-
- def estimate_apvoigt(self, x, y):
- """Estimation of *Area, Position, FWHM1, eta* of peaks, for
- pseudo-Voigt curves.
-
- This functions uses :meth:`estimate_pvoigt`, then converts the height
- parameter to area.
-
- :param x: Array of abscissa values
- :param y: Array of ordinate values (``y = f(x)``)
- :return: Tuple of estimated fit parameters and fit constraints.
- Parameters to be estimated for each peak are:
- *Area, Position, FWHM, eta*.
- """
- fittedpar, cons = self.estimate_pvoigt(x, y)
- npeaks = len(fittedpar) // 4
- # Assume 50% of the area is determined by the gaussian and 50% by
- # the Lorentzian.
- for i in range(npeaks):
- height = fittedpar[4 * i]
- fwhm = fittedpar[4 * i + 2]
- fittedpar[4 * i] = 0.5 * (height * fwhm * 0.5 * numpy.pi) +\
- 0.5 * (height * fwhm / (2.0 * numpy.sqrt(2 * numpy.log(2)))
- ) * numpy.sqrt(2 * numpy.pi)
- return fittedpar, cons
-
- def estimate_ahypermet(self, x, y):
- """Estimation of *area, position, fwhm, st_area_r, st_slope_r,
- lt_area_r, lt_slope_r, step_height_r* of peaks, for hypermet curves.
-
- :param x: Array of abscissa values
- :param y: Array of ordinate values (``y = f(x)``)
- :return: Tuple of estimated fit parameters and fit constraints.
- Parameters to be estimated for each peak are:
- *area, position, fwhm, st_area_r, st_slope_r,
- lt_area_r, lt_slope_r, step_height_r* .
- """
- yscaling = self.config.get('Yscaling', 1.0)
- if yscaling == 0:
- yscaling = 1.0
- fittedpar, cons = self.estimate_height_position_fwhm(x, y)
- npeaks = len(fittedpar) // 3
- newpar = []
- newcons = numpy.zeros((8 * npeaks, 3), numpy.float64)
- main_peak = 0
- # find out related parameters proper index
- if not self.config['NoConstraintsFlag']:
- if self.config['SameFwhmFlag']:
- j = 0
- # get the index of the free FWHM
- for i in range(npeaks):
- if cons[3 * i + 2, 0] != 4:
- j = i
- for i in range(npeaks):
- if i != j:
- cons[3 * i + 2, 1] = 8 * j + 2
- main_peak = j
- for i in range(npeaks):
- if fittedpar[3 * i] > fittedpar[3 * main_peak]:
- main_peak = i
-
- for i in range(npeaks):
- height = fittedpar[3 * i]
- position = fittedpar[3 * i + 1]
- fwhm = fittedpar[3 * i + 2]
- area = (height * fwhm / (2.0 * numpy.sqrt(2 * numpy.log(2)))
- ) * numpy.sqrt(2 * numpy.pi)
- # the gaussian parameters
- newpar.append(area)
- newpar.append(position)
- newpar.append(fwhm)
- # print "area, pos , fwhm = ",area,position,fwhm
- # Avoid zero derivatives because of not calculating contribution
- g_term = 1
- st_term = 1
- lt_term = 1
- step_term = 1
- if self.config['HypermetTails'] != 0:
- g_term = self.config['HypermetTails'] & 1
- st_term = (self.config['HypermetTails'] >> 1) & 1
- lt_term = (self.config['HypermetTails'] >> 2) & 1
- step_term = (self.config['HypermetTails'] >> 3) & 1
- if g_term == 0:
- # fix the gaussian parameters
- newcons[8 * i, 0] = CFIXED
- newcons[8 * i + 1, 0] = CFIXED
- newcons[8 * i + 2, 0] = CFIXED
- # the short tail parameters
- if ((area * yscaling) <
- self.config['MinGaussArea4ShortTail']) | \
- (st_term == 0):
- newpar.append(0.0)
- newpar.append(0.0)
- newcons[8 * i + 3, 0] = CFIXED
- newcons[8 * i + 3, 1] = 0.0
- newcons[8 * i + 3, 2] = 0.0
- newcons[8 * i + 4, 0] = CFIXED
- newcons[8 * i + 4, 1] = 0.0
- newcons[8 * i + 4, 2] = 0.0
- else:
- newpar.append(self.config['InitialShortTailAreaRatio'])
- newpar.append(self.config['InitialShortTailSlopeRatio'])
- newcons[8 * i + 3, 0] = CQUOTED
- newcons[8 * i + 3, 1] = self.config['MinShortTailAreaRatio']
- newcons[8 * i + 3, 2] = self.config['MaxShortTailAreaRatio']
- newcons[8 * i + 4, 0] = CQUOTED
- newcons[8 * i + 4, 1] = self.config['MinShortTailSlopeRatio']
- newcons[8 * i + 4, 2] = self.config['MaxShortTailSlopeRatio']
- # the long tail parameters
- if ((area * yscaling) <
- self.config['MinGaussArea4LongTail']) | \
- (lt_term == 0):
- newpar.append(0.0)
- newpar.append(0.0)
- newcons[8 * i + 5, 0] = CFIXED
- newcons[8 * i + 5, 1] = 0.0
- newcons[8 * i + 5, 2] = 0.0
- newcons[8 * i + 6, 0] = CFIXED
- newcons[8 * i + 6, 1] = 0.0
- newcons[8 * i + 6, 2] = 0.0
- else:
- newpar.append(self.config['InitialLongTailAreaRatio'])
- newpar.append(self.config['InitialLongTailSlopeRatio'])
- newcons[8 * i + 5, 0] = CQUOTED
- newcons[8 * i + 5, 1] = self.config['MinLongTailAreaRatio']
- newcons[8 * i + 5, 2] = self.config['MaxLongTailAreaRatio']
- newcons[8 * i + 6, 0] = CQUOTED
- newcons[8 * i + 6, 1] = self.config['MinLongTailSlopeRatio']
- newcons[8 * i + 6, 2] = self.config['MaxLongTailSlopeRatio']
- # the step parameters
- if ((height * yscaling) <
- self.config['MinGaussHeight4StepTail']) | \
- (step_term == 0):
- newpar.append(0.0)
- newcons[8 * i + 7, 0] = CFIXED
- newcons[8 * i + 7, 1] = 0.0
- newcons[8 * i + 7, 2] = 0.0
- else:
- newpar.append(self.config['InitialStepTailHeightRatio'])
- newcons[8 * i + 7, 0] = CQUOTED
- newcons[8 * i + 7, 1] = self.config['MinStepTailHeightRatio']
- newcons[8 * i + 7, 2] = self.config['MaxStepTailHeightRatio']
- # if self.config['NoConstraintsFlag'] == 1:
- # newcons=numpy.zeros((8*npeaks, 3),numpy.float64)
- if npeaks > 0:
- if g_term:
- if self.config['PositiveHeightAreaFlag']:
- for i in range(npeaks):
- newcons[8 * i, 0] = CPOSITIVE
- if self.config['PositiveFwhmFlag']:
- for i in range(npeaks):
- newcons[8 * i + 2, 0] = CPOSITIVE
- if self.config['SameFwhmFlag']:
- for i in range(npeaks):
- if i != main_peak:
- newcons[8 * i + 2, 0] = CFACTOR
- newcons[8 * i + 2, 1] = 8 * main_peak + 2
- newcons[8 * i + 2, 2] = 1.0
- if self.config['HypermetQuotedPositionFlag']:
- for i in range(npeaks):
- delta = self.config['DeltaPositionFwhmUnits'] * fwhm
- newcons[8 * i + 1, 0] = CQUOTED
- newcons[8 * i + 1, 1] = newpar[8 * i + 1] - delta
- newcons[8 * i + 1, 2] = newpar[8 * i + 1] + delta
- if self.config['SameSlopeRatioFlag']:
- for i in range(npeaks):
- if i != main_peak:
- newcons[8 * i + 4, 0] = CFACTOR
- newcons[8 * i + 4, 1] = 8 * main_peak + 4
- newcons[8 * i + 4, 2] = 1.0
- newcons[8 * i + 6, 0] = CFACTOR
- newcons[8 * i + 6, 1] = 8 * main_peak + 6
- newcons[8 * i + 6, 2] = 1.0
- if self.config['SameAreaRatioFlag']:
- for i in range(npeaks):
- if i != main_peak:
- newcons[8 * i + 3, 0] = CFACTOR
- newcons[8 * i + 3, 1] = 8 * main_peak + 3
- newcons[8 * i + 3, 2] = 1.0
- newcons[8 * i + 5, 0] = CFACTOR
- newcons[8 * i + 5, 1] = 8 * main_peak + 5
- newcons[8 * i + 5, 2] = 1.0
- return newpar, newcons
-
- def estimate_stepdown(self, x, y):
- """Estimation of parameters for stepdown curves.
-
- The functions estimates gaussian parameters for the derivative of
- the data, takes the largest gaussian peak and uses its estimated
- parameters to define the center of the step and its fwhm. The
- estimated amplitude returned is simply ``max(y) - min(y)``.
-
- :param x: Array of abscissa values
- :param y: Array of ordinate values (``y = f(x)``)
- :return: Tuple of estimated fit parameters and fit newconstraints.
- Parameters to be estimated for each stepdown are:
- *height, centroid, fwhm* .
- """
- crappyfilter = [-0.25, -0.75, 0.0, 0.75, 0.25]
- cutoff = len(crappyfilter) // 2
- y_deriv = numpy.convolve(y,
- crappyfilter,
- mode="valid")
-
- # make the derivative's peak have the same amplitude as the step
- if max(y_deriv) > 0:
- y_deriv = y_deriv * max(y) / max(y_deriv)
-
- fittedpar, newcons = self.estimate_height_position_fwhm(
- x[cutoff:-cutoff], y_deriv)
-
- data_amplitude = max(y) - min(y)
-
- # use parameters from largest gaussian found
- if len(fittedpar):
- npeaks = len(fittedpar) // 3
- largest_index = 0
- largest = [data_amplitude,
- fittedpar[3 * largest_index + 1],
- fittedpar[3 * largest_index + 2]]
- for i in range(npeaks):
- if fittedpar[3 * i] > largest[0]:
- largest_index = i
- largest = [data_amplitude,
- fittedpar[3 * largest_index + 1],
- fittedpar[3 * largest_index + 2]]
- else:
- # no peak was found
- largest = [data_amplitude, # height
- x[len(x)//2], # center: middle of x range
- self.config["FwhmPoints"] * (x[1] - x[0])] # fwhm: default value
-
- # Setup constrains
- newcons = numpy.zeros((3, 3), numpy.float64)
- if not self.config['NoConstraintsFlag']:
- # Setup height constrains
- if self.config['PositiveHeightAreaFlag']:
- newcons[0, 0] = CPOSITIVE
- newcons[0, 1] = 0
- newcons[0, 2] = 0
-
- # Setup position constrains
- if self.config['QuotedPositionFlag']:
- newcons[1, 0] = CQUOTED
- newcons[1, 1] = min(x)
- newcons[1, 2] = max(x)
-
- # Setup positive FWHM constrains
- if self.config['PositiveFwhmFlag']:
- newcons[2, 0] = CPOSITIVE
- newcons[2, 1] = 0
- newcons[2, 2] = 0
-
- return largest, newcons
-
- def estimate_slit(self, x, y):
- """Estimation of parameters for slit curves.
-
- The functions estimates stepup and stepdown parameters for the largest
- steps, and uses them for calculating the center (middle between stepup
- and stepdown), the height (maximum amplitude in data), the fwhm
- (distance between the up- and down-step centers) and the beamfwhm
- (average of FWHM for up- and down-step).
-
- :param x: Array of abscissa values
- :param y: Array of ordinate values (``y = f(x)``)
- :return: Tuple of estimated fit parameters and fit constraints.
- Parameters to be estimated for each slit are:
- *height, position, fwhm, beamfwhm* .
- """
- largestup, cons = self.estimate_stepup(x, y)
- largestdown, cons = self.estimate_stepdown(x, y)
- fwhm = numpy.fabs(largestdown[1] - largestup[1])
- beamfwhm = 0.5 * (largestup[2] + largestdown[1])
- beamfwhm = min(beamfwhm, fwhm / 10.0)
- beamfwhm = max(beamfwhm, (max(x) - min(x)) * 3.0 / len(x))
-
- y_minus_bg = y - self.strip_bg(y)
- height = max(y_minus_bg)
-
- i1 = numpy.nonzero(y_minus_bg >= 0.5 * height)[0]
- xx = numpy.take(x, i1)
- position = (xx[0] + xx[-1]) / 2.0
- fwhm = xx[-1] - xx[0]
- largest = [height, position, fwhm, beamfwhm]
- cons = numpy.zeros((4, 3), numpy.float64)
- # Setup constrains
- if not self.config['NoConstraintsFlag']:
- # Setup height constrains
- if self.config['PositiveHeightAreaFlag']:
- cons[0, 0] = CPOSITIVE
- cons[0, 1] = 0
- cons[0, 2] = 0
-
- # Setup position constrains
- if self.config['QuotedPositionFlag']:
- cons[1, 0] = CQUOTED
- cons[1, 1] = min(x)
- cons[1, 2] = max(x)
-
- # Setup positive FWHM constrains
- if self.config['PositiveFwhmFlag']:
- cons[2, 0] = CPOSITIVE
- cons[2, 1] = 0
- cons[2, 2] = 0
-
- # Setup positive FWHM constrains
- if self.config['PositiveFwhmFlag']:
- cons[3, 0] = CPOSITIVE
- cons[3, 1] = 0
- cons[3, 2] = 0
- return largest, cons
-
- def estimate_stepup(self, x, y):
- """Estimation of parameters for a single step up curve.
-
- The functions estimates gaussian parameters for the derivative of
- the data, takes the largest gaussian peak and uses its estimated
- parameters to define the center of the step and its fwhm. The
- estimated amplitude returned is simply ``max(y) - min(y)``.
-
- :param x: Array of abscissa values
- :param y: Array of ordinate values (``y = f(x)``)
- :return: Tuple of estimated fit parameters and fit constraints.
- Parameters to be estimated for each stepup are:
- *height, centroid, fwhm* .
- """
- crappyfilter = [0.25, 0.75, 0.0, -0.75, -0.25]
- cutoff = len(crappyfilter) // 2
- y_deriv = numpy.convolve(y, crappyfilter, mode="valid")
- if max(y_deriv) > 0:
- y_deriv = y_deriv * max(y) / max(y_deriv)
-
- fittedpar, cons = self.estimate_height_position_fwhm(
- x[cutoff:-cutoff], y_deriv)
-
- # for height, use the data amplitude after removing the background
- data_amplitude = max(y) - min(y)
-
- # find params of the largest gaussian found
- if len(fittedpar):
- npeaks = len(fittedpar) // 3
- largest_index = 0
- largest = [data_amplitude,
- fittedpar[3 * largest_index + 1],
- fittedpar[3 * largest_index + 2]]
- for i in range(npeaks):
- if fittedpar[3 * i] > largest[0]:
- largest_index = i
- largest = [fittedpar[3 * largest_index],
- fittedpar[3 * largest_index + 1],
- fittedpar[3 * largest_index + 2]]
- else:
- # no peak was found
- largest = [data_amplitude, # height
- x[len(x)//2], # center: middle of x range
- self.config["FwhmPoints"] * (x[1] - x[0])] # fwhm: default value
-
- newcons = numpy.zeros((3, 3), numpy.float64)
- # Setup constrains
- if not self.config['NoConstraintsFlag']:
- # Setup height constraints
- if self.config['PositiveHeightAreaFlag']:
- newcons[0, 0] = CPOSITIVE
- newcons[0, 1] = 0
- newcons[0, 2] = 0
-
- # Setup position constraints
- if self.config['QuotedPositionFlag']:
- newcons[1, 0] = CQUOTED
- newcons[1, 1] = min(x)
- newcons[1, 2] = max(x)
-
- # Setup positive FWHM constraints
- if self.config['PositiveFwhmFlag']:
- newcons[2, 0] = CPOSITIVE
- newcons[2, 1] = 0
- newcons[2, 2] = 0
-
- return largest, newcons
-
- def estimate_periodic_gauss(self, x, y):
- """Estimation of parameters for periodic gaussian curves:
- *number of peaks, distance between peaks, height, position of the
- first peak, fwhm*
-
- The functions detects all peaks, then computes the parameters the
- following way:
-
- - *distance*: average of distances between detected peaks
- - *height*: average height of detected peaks
- - *fwhm*: fwhm of the highest peak (in number of samples) if
- field ``'AutoFwhm'`` in :attr:`config` is ``True``, else take
- the default value (field ``'FwhmPoints'`` in :attr:`config`)
-
- :param x: Array of abscissa values
- :param y: Array of ordinate values (``y = f(x)``)
- :return: Tuple of estimated fit parameters and fit constraints.
- """
- yscaling = self.config.get('Yscaling', 1.0)
- if yscaling == 0:
- yscaling = 1.0
-
- bg = self.strip_bg(y)
-
- if self.config['AutoFwhm']:
- search_fwhm = guess_fwhm(y)
- else:
- search_fwhm = int(float(self.config['FwhmPoints']))
- search_sens = float(self.config['Sensitivity'])
-
- if search_fwhm < 3:
- search_fwhm = 3
-
- if search_sens < 1:
- search_sens = 1
-
- if len(y) > 1.5 * search_fwhm:
- peaks = peak_search(yscaling * y, fwhm=search_fwhm,
- sensitivity=search_sens)
- else:
- peaks = []
- npeaks = len(peaks)
- if not npeaks:
- fittedpar = []
- cons = numpy.zeros((len(fittedpar), 3), numpy.float64)
- return fittedpar, cons
-
- fittedpar = [0.0, 0.0, 0.0, 0.0, 0.0]
-
- # The number of peaks
- fittedpar[0] = npeaks
-
- # The separation between peaks in x units
- delta = 0.0
- height = 0.0
- for i in range(npeaks):
- height += y[int(peaks[i])] - bg[int(peaks[i])]
- if i != npeaks - 1:
- delta += (x[int(peaks[i + 1])] - x[int(peaks[i])])
-
- # delta between peaks
- if npeaks > 1:
- fittedpar[1] = delta / (npeaks - 1)
-
- # starting height
- fittedpar[2] = height / npeaks
-
- # position of the first peak
- fittedpar[3] = x[int(peaks[0])]
-
- # Estimate the fwhm
- fittedpar[4] = search_fwhm
-
- # setup constraints
- cons = numpy.zeros((5, 3), numpy.float64)
- cons[0, 0] = CFIXED # the number of gaussians
- if npeaks == 1:
- cons[1, 0] = CFIXED # the delta between peaks
- else:
- cons[1, 0] = CFREE
- j = 2
- # Setup height area constrains
- if not self.config['NoConstraintsFlag']:
- if self.config['PositiveHeightAreaFlag']:
- # POSITIVE = 1
- cons[j, 0] = CPOSITIVE
- cons[j, 1] = 0
- cons[j, 2] = 0
- j += 1
-
- # Setup position constrains
- if not self.config['NoConstraintsFlag']:
- if self.config['QuotedPositionFlag']:
- # QUOTED = 2
- cons[j, 0] = CQUOTED
- cons[j, 1] = min(x)
- cons[j, 2] = max(x)
- j += 1
-
- # Setup positive FWHM constrains
- if not self.config['NoConstraintsFlag']:
- if self.config['PositiveFwhmFlag']:
- # POSITIVE=1
- cons[j, 0] = CPOSITIVE
- cons[j, 1] = 0
- cons[j, 2] = 0
- j += 1
- return fittedpar, cons
-
- def configure(self, **kw):
- """Add new / unknown keyword arguments to :attr:`config`,
- update entries in :attr:`config` if the parameter name is a existing
- key.
-
- :param kw: Dictionary of keyword arguments.
- :return: Configuration dictionary :attr:`config`
- """
- if not kw.keys():
- return self.config
- for key in kw.keys():
- notdone = 1
- # take care of lower / upper case problems ...
- for config_key in self.config.keys():
- if config_key.lower() == key.lower():
- self.config[config_key] = kw[key]
- notdone = 0
- if notdone:
- self.config[key] = kw[key]
- return self.config
-
-fitfuns = FitTheories()
-
-THEORY = OrderedDict((
- ('Gaussians',
- FitTheory(description='Gaussian functions',
- function=functions.sum_gauss,
- parameters=('Height', 'Position', 'FWHM'),
- estimate=fitfuns.estimate_height_position_fwhm,
- configure=fitfuns.configure)),
- ('Lorentz',
- FitTheory(description='Lorentzian functions',
- function=functions.sum_lorentz,
- parameters=('Height', 'Position', 'FWHM'),
- estimate=fitfuns.estimate_height_position_fwhm,
- configure=fitfuns.configure)),
- ('Area Gaussians',
- FitTheory(description='Gaussian functions (area)',
- function=functions.sum_agauss,
- parameters=('Area', 'Position', 'FWHM'),
- estimate=fitfuns.estimate_agauss,
- configure=fitfuns.configure)),
- ('Area Lorentz',
- FitTheory(description='Lorentzian functions (area)',
- function=functions.sum_alorentz,
- parameters=('Area', 'Position', 'FWHM'),
- estimate=fitfuns.estimate_alorentz,
- configure=fitfuns.configure)),
- ('Pseudo-Voigt Line',
- FitTheory(description='Pseudo-Voigt functions',
- function=functions.sum_pvoigt,
- parameters=('Height', 'Position', 'FWHM', 'Eta'),
- estimate=fitfuns.estimate_pvoigt,
- configure=fitfuns.configure)),
- ('Area Pseudo-Voigt',
- FitTheory(description='Pseudo-Voigt functions (area)',
- function=functions.sum_apvoigt,
- parameters=('Area', 'Position', 'FWHM', 'Eta'),
- estimate=fitfuns.estimate_apvoigt,
- configure=fitfuns.configure)),
- ('Split Gaussian',
- FitTheory(description='Asymmetric gaussian functions',
- function=functions.sum_splitgauss,
- parameters=('Height', 'Position', 'LowFWHM',
- 'HighFWHM'),
- estimate=fitfuns.estimate_splitgauss,
- configure=fitfuns.configure)),
- ('Split Lorentz',
- FitTheory(description='Asymmetric lorentzian functions',
- function=functions.sum_splitlorentz,
- parameters=('Height', 'Position', 'LowFWHM', 'HighFWHM'),
- estimate=fitfuns.estimate_splitgauss,
- configure=fitfuns.configure)),
- ('Split Pseudo-Voigt',
- FitTheory(description='Asymmetric pseudo-Voigt functions',
- function=functions.sum_splitpvoigt,
- parameters=('Height', 'Position', 'LowFWHM',
- 'HighFWHM', 'Eta'),
- estimate=fitfuns.estimate_splitpvoigt,
- configure=fitfuns.configure)),
- ('Step Down',
- FitTheory(description='Step down function',
- function=functions.sum_stepdown,
- parameters=('Height', 'Position', 'FWHM'),
- estimate=fitfuns.estimate_stepdown,
- configure=fitfuns.configure)),
- ('Step Up',
- FitTheory(description='Step up function',
- function=functions.sum_stepup,
- parameters=('Height', 'Position', 'FWHM'),
- estimate=fitfuns.estimate_stepup,
- configure=fitfuns.configure)),
- ('Slit',
- FitTheory(description='Slit function',
- function=functions.sum_slit,
- parameters=('Height', 'Position', 'FWHM', 'BeamFWHM'),
- estimate=fitfuns.estimate_slit,
- configure=fitfuns.configure)),
- ('Atan',
- FitTheory(description='Arctan step up function',
- function=functions.atan_stepup,
- parameters=('Height', 'Position', 'Width'),
- estimate=fitfuns.estimate_stepup,
- configure=fitfuns.configure)),
- ('Hypermet',
- FitTheory(description='Hypermet functions',
- function=fitfuns.ahypermet, # customized version of functions.sum_ahypermet
- parameters=('G_Area', 'Position', 'FWHM', 'ST_Area',
- 'ST_Slope', 'LT_Area', 'LT_Slope', 'Step_H'),
- estimate=fitfuns.estimate_ahypermet,
- configure=fitfuns.configure)),
- # ('Periodic Gaussians',
- # FitTheory(description='Periodic gaussian functions',
- # function=functions.periodic_gauss,
- # parameters=('N', 'Delta', 'Height', 'Position', 'FWHM'),
- # estimate=fitfuns.estimate_periodic_gauss,
- # configure=fitfuns.configure))
- ('Degree 2 Polynomial',
- FitTheory(description='Degree 2 polynomial'
- '\ny = a*x^2 + b*x +c',
- function=fitfuns.poly,
- parameters=['a', 'b', 'c'],
- estimate=fitfuns.estimate_quadratic)),
- ('Degree 3 Polynomial',
- FitTheory(description='Degree 3 polynomial'
- '\ny = a*x^3 + b*x^2 + c*x + d',
- function=fitfuns.poly,
- parameters=['a', 'b', 'c', 'd'],
- estimate=fitfuns.estimate_cubic)),
- ('Degree 4 Polynomial',
- FitTheory(description='Degree 4 polynomial'
- '\ny = a*x^4 + b*x^3 + c*x^2 + d*x + e',
- function=fitfuns.poly,
- parameters=['a', 'b', 'c', 'd', 'e'],
- estimate=fitfuns.estimate_quartic)),
- ('Degree 5 Polynomial',
- FitTheory(description='Degree 5 polynomial'
- '\ny = a*x^5 + b*x^4 + c*x^3 + d*x^2 + e*x + f',
- function=fitfuns.poly,
- parameters=['a', 'b', 'c', 'd', 'e', 'f'],
- estimate=fitfuns.estimate_quintic)),
-))
-"""Dictionary of fit theories: fit functions and their associated estimation
-function, parameters list, configuration function and description.
-"""
-
-
-def test(a):
- from silx.math.fit import fitmanager
- x = numpy.arange(1000).astype(numpy.float64)
- p = [1500, 100., 50.0,
- 1500, 700., 50.0]
- y_synthetic = functions.sum_gauss(x, *p) + 1
-
- fit = fitmanager.FitManager(x, y_synthetic)
- fit.addtheory('Gaussians', functions.sum_gauss, ['Height', 'Position', 'FWHM'],
- a.estimate_height_position_fwhm)
- fit.settheory('Gaussians')
- fit.setbackground('Linear')
-
- fit.estimate()
- fit.runfit()
-
- y_fit = fit.gendata()
-
- print("Fit parameter names: %s" % str(fit.get_names()))
- print("Theoretical parameters: %s" % str(numpy.append([1, 0], p)))
- print("Fitted parameters: %s" % str(fit.get_fitted_parameters()))
-
- try:
- from silx.gui import qt
- from silx.gui.plot import plot1D
- app = qt.QApplication([])
-
- # Offset of 1 to see the difference in log scale
- plot1D(x, (y_synthetic + 1, y_fit), "Input data + 1, Fit")
-
- app.exec_()
- except ImportError:
- _logger.warning("Unable to load qt binding, can't plot results.")
-
-
-if __name__ == "__main__":
- test(fitfuns)
generated by cgit v1.2.3 (git 2.25.1) at 2024-04-23 23:01:59 +0000