qdsfit/QDS.py

394 lines
15 KiB
Python
Executable File

#!/usr/bin/env python
# -*- encoding: utf-8 -*-
import os
import sys
import re
import signal
from PyQt4.QtCore import *
from PyQt4.QtGui import *
import matplotlib
from mathlib import fit_anneal, fit_lbfgsb, fit_odr_cmplx
matplotlib.use('agg')
from matplotlib import pyplot
from matplotlib.colors import hex2color
#matplotlib.rc_file("default.mplrc")
import numpy as N
import QDSMain
from data import Data, Conductivity, Peak
import pyqtgraph as pg
#import yaff
from CustomWidgets import ParameterWidget
USE_CROSSH=False
class AppWindow(QMainWindow):
def __init__(self, parent=None):
super(AppWindow, self).__init__(parent)
self.ui = QDSMain.Ui_MainWindow()
self.ui.setupUi(self)
self.picked_artist = None
self.data = None
self.Conductivity = None
self._lines = dict()
self.peakId = 0
self.peakBoxes = {}
self.parameterWidget = ParameterWidget()
self.ui.dockWidget_3.setWidget(self.parameterWidget)
## menubar
fileMenu = self.menuBar().addMenu("File")
openFile = QAction("&Open", self)
openFile.setShortcut(QKeySequence.Open)
openFile.triggered.connect(self.openFile)
fileMenu.addAction(openFile)
saveFile = QAction("&Save Fit Result", self)
saveFile.setShortcut(QKeySequence.Save)
saveFile.triggered.connect(self.saveFitResult)
fileMenu.addAction(saveFile)
# fitting methods
fitMenu = self.menuBar().addMenu("Standard Fits")
# lm
fit_lmAction = QAction("&Levenberg-Marquardt", self)
fit_lmAction.setShortcut(QKeySequence("Ctrl+F"))
fitMenu.addAction(fit_lmAction)
# lbfgsb
fit_lbfgsbAction = QAction("&L-BFGS-B", self)
fitMenu.addAction(fit_lbfgsbAction)
# Simulated Annealing
fit_annealAction = QAction("&Simulated Annealing", self)
fitMenu.addAction(fit_annealAction)
# YAFF
yaffMenu = self.menuBar().addMenu("YAFF")
start_yaff = QAction("&Startparam.", self)
yaffMenu.addAction(start_yaff)
fit_yaff = QAction("&Fit", self)
yaffMenu.addAction(fit_yaff)
self.signalMapper = QSignalMapper(self)
for i, fit_action in enumerate([fit_lmAction, fit_lbfgsbAction, fit_annealAction
]):
self.signalMapper.setMapping(fit_action, i)
fit_action.triggered.connect(self.signalMapper.map)
self.signalMapper.mapped.connect(self.fitData)
# save fitted values
#self.ui.actionSave_FitResult.triggered.connect(self.saveFitResult)# replaced by menu
self.data = Data()
self.fit_boundary = pg.LinearRegionItem(brush=QColor(254,254,254,10))
self.ui.graphicsView.addItem(self.data.data_curve_imag)
self.ui.graphicsView.addItem(self.data.data_curve_real)
self.ui.graphicsView.addItem(self.data.fitted_curve_imag)
self.ui.graphicsView.addItem(self.data.fitted_curve_real)
# cross hair
if USE_CROSSH:
#self.poslabel = pg.LabelItem(justify='right')
self.horizontalLine = pg.InfiniteLine(angle=0)
self.verticalLine = pg.InfiniteLine(angle=90)
#self.ui.graphicsView.plotItem.addItem(self.poslabel)
self.ui.graphicsView.addItem(self.horizontalLine)
self.ui.graphicsView.addItem(self.verticalLine)
self.ui.graphicsView.addItem(self.fit_boundary)
self.ui.graphicsView.setLogMode(x=True, y=True)
self.ui.graphicsView.showGrid(x=True, y=True)
self.ui.graphicsView.setLabel("bottom","Frequency",units="Hz")
self.ui.graphicsView.setLabel("left", "Dielectric loss", units="Debye")
self.ui.graphicsView.disableAutoRange()
#self.peak_box_layout = QGridLayout()
#self.ui.scrollAreaWidgetContents.setLayout(self.peak_box_layout)
##self.huh = pg.SignalProxy(, slot=self.mPE)
sc = self.ui.graphicsView.scene()
sc.sigMouseClicked.connect(self.mousePress)
sc.sigMouseMoved.connect(self.updateCrosshair)
def updateCrosshair(self,evt):
vb = self.ui.graphicsView.getPlotItem().vb
pos = vb.mapSceneToView(evt)
self.last_pos = pos
if USE_CROSSH:
self.horizontalLine.setBounds([self.data.frequency.min(), self.data.frequency.max()])
self.horizontalLine.setPos(pos.y())
self.verticalLine.setPos(pos.x())
#self.poslabel.setText("<span style='font-size: 12pt'>f=%0.1f Hz <span style='color: red'>e=%0.1f</span>"\
# % (pos.x(),pos.y()))
def mousePress(self, evt):
vb = self.ui.graphicsView.getPlotItem().vb
data_pos = self.last_pos
if self.ui.actionAdd_Peak.isChecked() and not self.ui.actionAdd_Cond.isChecked():
self.addPeak(data_pos)
self.ui.actionAdd_Peak.setChecked(False)
if self.ui.actionAdd_Cond.isChecked() and not self.ui.actionAdd_Peak.isChecked():
self.addCond(data_pos)
self.ui.actionAdd_Cond.setChecked(False)
def saveFitResult(self):
"""
Saving fit parameters to fitresults.log
including temperature
"""
self.saveFitFigure()
if not os.path.exists("fitresults.log"):
f = open("fitresults.log", "w")
else:
f = open("fitresults.log", "a")
# write header
f.write("#%7s"%('T'))
parfmt = "%8.2f" # T formatting
# if self.Conductivity != None: pass# always true
f.write("%9s%9s%9s " % ("e_s", "sig", "pow_sig"))
parfmt += "%9.3g%9.3g%9.2f " # conductivity formatting
for i, pb in enumerate(self.peakBoxes):
enum_peak = ("e_inf_%i" % i, "tau_%i" % i, "alpha_%i" % i, "beta_%i" % i)
f.write("%9s%9s%9s%9s " % enum_peak)
print enum_peak
parfmt += "%9.3g%9.3g%9.2f%9.2f" # peak formatting
f.write("fit_xlow fit_xhigh") # TODO: store limits
parfmt += "%9.3g%9.3g"
f.write('\n')
f.flush()
#f.write("%3.2f "%(self.data.meta["T"]))
pars = list(self.fitresult)
pars.insert(0, self.data.meta["T"])
pars.append(self.data.fit_limits[0])
pars.append(self.data.fit_limits[1])
N.savetxt(f, N.array([pars, ]), fmt=parfmt, delimiter=" ")
f.close()
def saveFitFigure(self):
fig = pyplot.figure(figsize=(3.54, 2.75))
font = {'family' : 'sans serif',
'weight' : 'normal',
'size' : 8}
matplotlib.rc('font', **font)
pyplot.loglog(self.data.frequency, self.data.epsilon.imag, 'bo', markersize=3, label="Data")
pyplot.loglog(self.data.frequency, self.data.epsilon_fit, 'r-', lw=1, label="Fit")
for i,peak in enumerate(self.peakBoxes):
f,eps = peak.get_data()
color = hex2color(str(peak.get_color().name()))
pyplot.loglog(f,eps, ls="--", color=color , lw=0.75, label="Peak %i"%i)
if self.Conductivity != None:
f,eps = self.Conductivity.get_conductivity()
color = hex2color(str(self.Conductivity.get_color().name()))
pyplot.loglog(f,eps, ls="-.", color=color, lw=0.75, label="Cond.")
f,eps = self.Conductivity.get_epsilon_static()
pyplot.loglog(f,eps, ls=":", color=color, lw=0.75, label=r'$\epsilon_0$')
for i in (0,1): pyplot.axvline(x=self.data.fit_limits[i], color='g', ls="--")
pyplot.legend(title = "T=%.1f K"%(self.data.meta["T"]))
pyplot.grid()
pyplot.xlabel('f/Hz')
pyplot.ylabel('eps"')
pyplot.savefig(os.path.splitext(self.filepath)[0]+".png")
fig.clear()
def addCond(self, pos):
win = self.parameterWidget
win.setWindowFlags(win.windowFlags() | Qt.WindowStaysOnTopHint)
if self.Conductivity != None:
return
self.statusBar().showMessage("Click on graph")
self.Conductivity = Conductivity(mpl=self.ui.graphicsView, limits=self.data.fit_limits)
self.Conductivity.changedData.connect(self.updatePlot)
#print pos.x(), pos.y(),"test", 10**pos.x(),10**pos.y()
cond_par = [0, 10**(pos.y() + pos.x())*2*N.pi , 1.0]
self.Conductivity.setParameter(beta=cond_par)
self.parameterWidget.vlayout.insertWidget(1,self.Conductivity.widget)
self.Conductivity.widget.ui.removeButton.clicked.connect(self.delCond)
self.parameterWidget.vlayout.update()
self.parameterWidget.showNormal()
self.parameterWidget.raise_()
def delCond(self):
self.cond_param = None
self.cond = None
self.ui.graphicsView.removeItem(self.Conductivity.mpl_line)
self.ui.graphicsView.removeItem(self.Conductivity.mpl_line_static)
del self.Conductivity
self.Conductivity = None
self.updatePlot()
def addPeak(self, pos):
win = self.parameterWidget
win.setWindowFlags(win.windowFlags() | Qt.WindowStaysOnTopHint)
id_list = [ self.peakBoxes[key] for key in self.peakBoxes.keys()]
self.peakId = 1
while self.peakId in id_list:
self.peakId += 1
self.statusBar().showMessage("Select Peak Position")
peak = Peak(id=self.peakId, mpl=self.ui.graphicsView, limits=self.data.fit_limits)
# connect to delPeak
peak.widget.ui.removeButton.clicked.connect(self.delPeak)
#peak.widget.ui.removeButton.clicked.connect(self.delPeak)
peak.changedData.connect(self.updatePlot)
# delta_eps=(10**pos.y())*2, tau=1 / (2*N.pi*10**pos.x()), a=1, b=1
new_peak = [2*10**pos.y(), 1 / (2*N.pi*10**pos.x()), 1, 1]
peak.setParameter(beta = new_peak)
self.peakBoxes[peak] = self.peakId
for i,pb in enumerate(self.peakBoxes.keys()):
self.parameterWidget.vlayout.insertWidget(1,pb.widget)
self.parameterWidget.vlayout.update()
self.parameterWidget.showNormal()
self.parameterWidget.raise_()
def delPeak(self):
for i, peak in enumerate(self.peakBoxes.keys()):
if peak.widget.isHidden():
self.ui.graphicsView.removeItem(peak.mpl_line)
self.parameterWidget.vlayout.removeWidget(peak.widget)
self.peakBoxes.pop(peak)
self.parameterWidget.vlayout.update()
self.updatePlot()
def fitData(self, method):
if self.Conductivity != None:
start_parameter = list(self.Conductivity.getParameter())
fixed_params = [i for i in self.Conductivity.getFixed()]
else:
start_parameter = [0, 0, 1]
fixed_params = [0, 0, 0]
for pb in self.peakBoxes.keys():
[start_parameter.append(i) for i in pb.getParameter()]
[fixed_params.append(i) for i in pb.getFixed()]
log10fmin, log10fmax = self.fit_boundary.getRegion()
self.data.set_fit_xlimits(10**log10fmin, 10**log10fmax)
fit_methods = [fit_odr_cmplx, fit_lbfgsb, fit_anneal]
print "StartParameter", start_parameter
print "FixedParameter", fixed_params
print "Limits (xmin, xmax, ymin, ymax)", self.data.fit_limits
_freq, _fit = self.data.get_data()
# build function list
funcs = ["static","power"] if self.Conductivity != None else []
for pb in self.peakBoxes.keys():
funcs.append("hn")
newres = fit_methods[method](_freq, _fit, start_parameter, fixed_params, funcs)
print "Set fit data"
self.data.set_fit(newres.beta, funcs)
self.ui.statusbar.showMessage(" ".join(newres.stopreason))
result = newres.beta
self.fitresult = result
for i, pb in enumerate(self.peakBoxes.keys()):
beta_peak = newres.beta [3 + i*4 : 3 + (i + 1)*4]
sd_beta_peak = newres.sd_beta[3 + i*4 : 3 + (i + 1)*4]
pb.setParameter(beta = beta_peak, sd_beta = sd_beta_peak)
if self.Conductivity != None:
self.Conductivity.setParameter(beta = newres.beta[:3], sd_beta = newres.sd_beta[:3])
#print "*** FIT RESULTS ***"
#print u"\u03c3"
#print u"\u0394\u03b5"
self.updatePlot()
def openFile(self):
#path = unicode(QFileDialog.getOpenFileName(self, "Open file"))
path = "MCM42PG0_199.96K.dat"
self.filepath=path
# TODO analyize file (LF,MF, HF) and act accordingly
data = N.loadtxt(path, skiprows=4)
self.setWindowTitle(os.path.basename(path))
numpat = re.compile('\d+\.\d+')
try:
Temp = None
for line in open(path).readlines():
if re.search("Fixed", line) or re.search("Temp", line):
print "Found line with Fixed or Temp"
Temp = float(re.search(numpat, line).group())
print "Temperature found in file:", Temp
break
if Temp == None: raise ValueError
except:
Temp = QInputDialog.getDouble(self, "No temperature found in data set", "Temperature/K:", value=0.00)[0]
# mask the data to values > 0 (loglog plot)
mask = (data[:, 1] > 0) & (data[:, 2] > 0) #& (data[:,2]>1e-3) & (data[:,0] > 1e-2)
_freq = data[mask, 0]
_die_stor = data[mask, 1]
_die_loss = data[mask, 2]
self.data.set_data(_freq, _die_stor, _die_loss)
self.data.meta["T"] = Temp
self.fit_boundary.setRegion([N.log10(_freq.min()), N.log10(_freq.max())])
self.ui.graphicsView.enableAutoRange()
self.updatePlot()
self.ui.graphicsView.disableAutoRange()
def updatePlot(self):
# assemble fit function
funcs = ["static", "power"]
p0 = self.Conductivity.getParameter() if self.Conductivity is not None else [0.0, 0.0, 1.0]
for peak in self.peakBoxes.keys():
params = peak.getParameter()
p0.extend(params)
funcs.append("hn")
# calculate fit
self.data.set_fit(p0, funcs)
# replot data and fit, TODO: replot only if data changed
self.data.data_curve_real.setData(self.data.frequency, self.data.epsilon.real)
self.data.data_curve_imag.setData(self.data.frequency, self.data.epsilon.imag)
if len(self.peakBoxes) > 0 or self.Conductivity != None:
self.data.fitted_curve_real.setData(self.data.frequency, self.data.epsilon_fit.real)
self.data.fitted_curve_imag.setData(self.data.frequency, self.data.epsilon_fit.imag)
def sigint_handler(*args):
"""Handler for the SIGINT signal (CTRL + C).
"""
sys.stderr.write('\r')
if QMessageBox.question(None, '', "Are you sure you want to quit?",
QMessageBox.Yes | QMessageBox.No,
QMessageBox.Yes) == QMessageBox.Yes:
QApplication.quit()
if __name__ == '__main__':
signal.signal(signal.SIGINT, sigint_handler)
app = QApplication(sys.argv)
timer = QTimer()
timer.start(1000) # Check every second for Strg-c on Cmd line
timer.timeout.connect(lambda: None)
main = AppWindow()
main.showMaximized()
main.raise_()
sys.exit(app.exec_())