Source code for pyspex.optimize
#!/usr/bin/env python
from . import pyspex_f2py
import math
[docs]class Fit:
"""This is the parent class of all fit related commands.
:ivar stat: Fit statistics.
:vartype stat: basestring
:ivar cstat: C-statistics value.
:vartype cstat: float
:ivar chisq: Chi-squared value.
:vartype chisq: float
:ivar wstat: W-statistics value (not recommended).
:vartype wstat: float
:ivar nfree: Degrees of freedom.
:vartype nfree: int
:ivar cstatexp: Expected C-statistics value (C-stat only).
:vartype cstatexp: float
:ivar cstatrms: RMS uncertainty on expected C-stat value.
:vartype cstatrms: float
:ivar ann_rt: Simulated annealing: rt
:vartype ann_rt: float
:ivar ann_eps: Simulated annealing: eps
:vartype ann_eps: float
:ivar ann_t: Simulated annealing: t
:vartype ann_t: float
:ivar ann_vm: Simulated annealing: vm
:vartype ann_vm: float
:ivar ann_ns: Simulated annealing: ns
:vartype ann_ns: int
:ivar ann_max: Simulated annealing: max evaluations
:vartype ann_max: int
:ivar ann_print: Simulated annealing: print flag
:vartype ann_print: int
"""
def __init__(self):
self.stat = 'csta' # Fit statistics
self.meth = 'clas' # Fit method
self.cstat = 0. # C-statistics value
self.chisq = 0. # Chi-square value
self.wstat = 0. # W-statistics value
self.nfree = 0 # Degrees of freedom
self.cstatexp = 0. # Expected C-statistics value (C-stat only)
self.cstatrms = 0. # RMS uncertainty on expected C-stat value
self.ann_rt = 0. # Simulated annealing: rt
self.ann_eps = 0. # Simulated annealing: eps
self.ann_t = 0. # Simulated annealing: t
self.ann_vm = 0. # Simulated annealing: vm
self.ann_ns = 0 # Simulated annealing: ns
self.ann_max = 0 # Simulated annealing: max evaluations
self.ann_print = 0 # Print flag
# Set 'fit print 1' automatically
self.print(True)
[docs] def fit(self, niter=100):
"""Execute the SPEX fit command. The maximum number of iterations (niter) can be optionally set.
:param niter: Number of fit iterations.
:type niter: int
"""
pyspex_f2py.api_optimize_fit.api_fit(iter=float(niter))
self.update()
[docs] def set_statistic(self, stat):
"""Set the desired fit statistics.
:param stat: Abbreviation for the fit statistics to be used. For example: 'csta', 'chi2', 'wsta'.
:type stat: str
"""
pyspex_f2py.api_optimize_fit.api_fit_set_statistics(stat)
[docs] def set_statistic_inst(self, stat, inst, reg):
"""Set the desired fit statistics per instrument and region.
:param stat: Abbreviation for the fit statistics to be used. For example: 'csta', 'chi2', 'wsta'.
:type stat: str
:param inst: Instrument number.
:type inst: int
:param reg: Region number.
:type reg: int
"""
pyspex_f2py.api_optimize_fit.api_fit_set_statistics_inst(stat, inst, reg)
[docs] def get_statistic(self):
"""Get the current type of statistics being used in the fit, for example chi2, cstat or wstat."""
self.stat = pyspex_f2py.api_optimize_fit.fit_stat.view('S4').tobytes().decode('utf-8')
[docs] def set_method(self, meth):
"""Set the desired fit statistics.
:param stat: Abbreviation for the fit statistics to be used. For example: 'csta', 'chi2', 'wsta'.
:type stat: str
"""
pyspex_f2py.api_optimize_fit.api_fit_set_method(meth)
[docs] def get_method(self):
"""Get the current type of statistics being used in the fit, for example chi2, cstat or wstat."""
self.meth = pyspex_f2py.api_optimize_fit.fit_meth.view('S4').tobytes().decode('utf-8')
[docs] def set_ann(self, param, value):
"""Set the simulating annealing method parameters.
:param param: Type of annealing parameter (rt, t, eps, vm, ns, max, or print).
:type param: str
:param value: Value of the parameter (will be converted to the nearest int if necessary).
:type value: float
"""
pyspex_f2py.api_optimize_fit.api_fit_set_ann(param, value)
self.update()
[docs] def print(self, status):
"""Print each fit iteration to the console (default is True). The status variable can be either True or False,
which means printing is on or off, respectively.
:param status: Set fit to high verbosity (True is yes, False is no).
:type status: bool
"""
if status:
pyspex_f2py.api_command('fit print 1')
else:
pyspex_f2py.api_command('fit print 0')
[docs] def update(self):
"""Get the most recent statistics values from SPEX."""
pyspex_f2py.api_optimize_fit.api_fit_update()
self.cstat = pyspex_f2py.api_optimize_fit.fit_cstat
self.chisq = pyspex_f2py.api_optimize_fit.fit_chisq
self.wstat = pyspex_f2py.api_optimize_fit.fit_wstat
self.nfree = pyspex_f2py.api_optimize_fit.fit_nfree
self.cstatexp = pyspex_f2py.api_optimize_fit.fit_cstatexp
self.cstatrms = pyspex_f2py.api_optimize_fit.fit_cstatrms
self.ann_rt = float(pyspex_f2py.api_optimize_fit.fit_ann_rt)
self.ann_eps = float(pyspex_f2py.api_optimize_fit.fit_ann_eps)
self.ann_t = float(pyspex_f2py.api_optimize_fit.fit_ann_t)
self.ann_vm = float(pyspex_f2py.api_optimize_fit.fit_ann_vm)
self.ann_ns = int(pyspex_f2py.api_optimize_fit.fit_ann_ns)
self.ann_max = int(pyspex_f2py.api_optimize_fit.fit_ann_maxevl)
self.ann_print = int(pyspex_f2py.api_optimize_fit.fit_ann_print)
self.get_method()
self.get_statistic()
[docs] def show(self):
"""Print the fit statistics to the terminal."""
print("---------------------------------------------------------------------------------------------------")
if self.meth == 'clas':
print("Fit method : Classical Levenberg-Marquardt")
elif self.meth == 'simp':
print("Fit method : Simplex")
elif self.meth == 'sian':
print("Fit method : Simulated annealing method")
else:
print("Fit method : {0:s}".format(self.meth))
if self.stat == 'csta':
print("Fit statistic : C-statistic")
elif self.stat == 'chi2':
print("Fit statistic : Chi-squared statistic")
elif self.stat == 'wsta':
print("Fit statistic : W-statistic")
else:
print("Fit statistic : {0:s}".format(self.stat))
print("C-statistic : {0:#7.2f}".format(self.cstat))
print("Expected C-stat : {0:#7.2f} +/- {1:#7.2f}".format(self.cstatexp, self.cstatrms))
print("Chi-squared value : {0:#7.2f}".format(self.chisq))
print("Degrees of freedom: {0:d}".format(self.nfree))
print("W-statistic : {0:#7.2f}".format(self.wstat))
[docs]class Error:
"""Class to calculate errors for free fit parameters.
:ivar sector: Sector number of parameter
:vartype sector: int
:ivar component: Component number of parameter
:vartype component: int
:ivar parameter: Parameter name
:vartype parameter: str
:ivar value: Parameter value
:vartype value: float
:ivar lerr: Lower error boundary
:vartype lerr: float
:ivar uerr: Upper error boundary
:vartype uerr: float
:ivar lc: Is there a lower C-stat or chi**2 value found?
:vartype lc: bool
:ivar cmin: Lowest C-stat or chi**2 value
:vartype cmin: float
:ivar pmin: Parameter value for which a better C-stat or Chi**2 was found
:vartype pmin: float
:ivar dchi: Delta C-stat or chi**2 to optimize for
:vartype dchi: float
:ivar calculated: Is the error calculated?
:vartype calculated: bool
"""
def __init__(self):
self.sector = 1 # Sector number of parameter
self.component = 1 # Component number of parameter
self.parameter = '' # Parameter name
self.value = 0.0 # Parameter value
self.lerr = 0.0 # Lower error boundary
self.uerr = 0.0 # Upper error boundary
self.lc = False # Is there a lower C-stat or chi**2 value found
self.cmin = 0.0 # Lowest C-stat or chi**2 value
self.pmin = 0.0 # Parameter value for which a better C-stat or Chi**2 was found
self.dchi = 1.0 # Delta C-stat or chi**2 to optimize for
self.calculated = False # Is the error calculated?
[docs] def error(self, isect, icomp, name, dchi=None):
"""Calculate the error value for a particular parameter.
:param isect: Sector number of the parameter.
:type isect: int
:param icomp: Component number of the parameter.
:type icomp: int
:param name: Parameter name.
:type name: str
:param dchi: (Optional) :math:`\Delta\chi^2` value to optimize for (Default: 1.0, 68% errors)
:type dchi: float
"""
if dchi is not None:
ier = self.set_dchi(dchi)
if ier != 0:
print("Error setting dchi for the error calculation.")
return ier
self.sector = isect
self.component = icomp
self.parameter = name
self.calculated = True
self.lc = False
pyspex_f2py.api_optimize_error.api_error_calc(float(isect), float(icomp), name)
self.value = float(pyspex_f2py.api_optimize_error.error_par)
self.lerr = float(pyspex_f2py.api_optimize_error.error_low)
self.uerr = float(pyspex_f2py.api_optimize_error.error_up)
if pyspex_f2py.api_optimize_error.error_lc:
self.lc = True
self.cmin = float(pyspex_f2py.api_optimize_error.error_cmin)
self.pmin = float(pyspex_f2py.api_optimize_error.error_pmin)
[docs] def set_dchi(self, dchi):
"""Set the delta c-stat or delta chi**2 value that the error search should optimize for.
The default value is 1.0.
:param dchi: :math:`\Delta\chi^2` value to optimize for (Default: 1.0, 68% errors)
:type dchi: float
"""
if (dchi > 0):
ier = pyspex_f2py.api_command('error dchi ' + str(dchi))
else:
print('Error: Delta C-stat or chi**2 needs to be positive.')
ier = 1
return ier
[docs] def get_value(self):
"""Convenience function to return the parameter value and the errors.
:return: A tuple with the parameter value, lower error and upper error.
:rtype: tuple
"""
if (self.calculated):
print('Returning parameter ' + str(self.sector) + ' ' + str(self.component) + ' ' + self.parameter)
return self.value, self.lerr, self.uerr
else:
print('Error: Error is not calculated yet.')
