source: trunk/testDeriv.py @ 1493

Last change on this file since 1493 was 1493, checked in by vondreele, 7 years ago

set move=True for GenAtom? call for FillUnitCell?
show bin width in calibration plot - also label items
add check if nothing refined in calibrate
fix crash at end of refinement for reflections with bad widths
fix testDeriv to show constrained variables
fix SCExtinction for constrained extinction coef

  • Property svn:eol-style set to native
  • Property svn:keywords set to Date Author Revision URL Id
File size: 7.8 KB
Line 
1# -*- coding: utf-8 -*-
2#testDeriv.py
3'''
4*testDeriv: Check derivative computation*
5=========================================
6
7Use this to check derivatives used in structure least squares
8refinement against numerical values computed in this script.
9
10To use set ``DEBUG=True`` in GSASIIstrMain.py (line 22, as of version
111110); run the least squares - one cycle is sufficient.  Do the "Save
12Results"; this will write the file testDeriv.dat in the local
13directory.
14
15Then run this program to see plots of derivatives for all
16parameters refined in the last least squares.  Shown will be numerical
17derivatives generated over all observations (including penalty terms)
18and the corresponding analytical ones produced in the least
19squares. They should match.
20'''
21
22import os
23import os.path as ospath
24import sys
25import time
26import cPickle
27import wx
28import numpy as np
29import matplotlib as mpl
30import GSASIIpath
31import GSASIIstrMath as G2stMth
32import GSASIItestplot as plot
33import GSASIImapvars as G2mv
34import pytexture as ptx
35ptx.pyqlmninit()            #initialize fortran arrays for spherical harmonics
36
37def create(parent):
38    return testDeriv(parent)
39   
40[wxID_FILEEXIT, wxID_FILEOPEN, wxID_MAKEPLOTS,
41] = [wx.NewId() for _init_coll_File_Items in range(3)]
42
43def FileDlgFixExt(dlg,file):            #this is needed to fix a problem in linux wx.FileDialog
44    ext = dlg.GetWildcard().split('|')[2*dlg.GetFilterIndex()+1].strip('*')
45    if ext not in file:
46        file += ext
47    return file
48   
49class testDeriv(wx.Frame):
50
51    def _init_ctrls(self, parent):
52        wx.Frame.__init__(self, name='testDeriv', parent=parent,
53            size=wx.Size(800, 250),style=wx.DEFAULT_FRAME_STYLE, title='Test Jacobian Derivatives')
54        self.testDerivMenu = wx.MenuBar()
55        self.File = wx.Menu(title='')
56        self.File.Append(help='Open testDeriv.dat', id=wxID_FILEOPEN,
57             kind=wx.ITEM_NORMAL,text='Open testDeriv.dat file')
58        self.File.Append(help='Make derivative plots',id=wxID_MAKEPLOTS,
59            kind=wx.ITEM_NORMAL,text='Make plots')
60        self.File.Append(help='Exit from testDeriv', id=wxID_FILEEXIT, kind=wx.ITEM_NORMAL,
61            text='Exit')
62        self.Bind(wx.EVT_MENU, self.OnTestRead, id=wxID_FILEOPEN)
63        self.Bind(wx.EVT_MENU,self.OnMakePlots,id=wxID_MAKEPLOTS)
64        self.Bind(wx.EVT_MENU, self.OnFileExit, id=wxID_FILEEXIT)
65        self.testDerivMenu.Append(menu=self.File, title='File')
66        self.SetMenuBar(self.testDerivMenu)
67        self.testDerivPanel = wx.ScrolledWindow(self)       
68        self.plotNB = plot.PlotNotebook()
69       
70    def __init__(self, parent):
71        self._init_ctrls(parent)
72        self.Bind(wx.EVT_CLOSE, self.ExitMain)   
73        self.dirname = ''
74        self.testfile = []
75        self.dataFrame = None
76
77    def ExitMain(self, event):
78        sys.exit()
79       
80    def OnFileExit(self,event):
81        if self.dataFrame:
82            self.dataFrame.Clear() 
83            self.dataFrame.Destroy()
84        self.Close()
85
86    def OnTestRead(self,event):
87        dlg = wx.FileDialog(self, 'Open testDeriv.dat file', '.', 'testDeriv.dat')
88        if self.dirname:
89            dlg.SetDirectory(self.dirname)
90        try:
91            if dlg.ShowModal() == wx.ID_OK:
92                self.dirname = dlg.GetDirectory()
93                testFile = dlg.GetPath()
94                file = open(testFile,'rb')
95                self.values = cPickle.load(file)
96                self.HistoPhases = cPickle.load(file)
97                (G2mv.dependentParmList,G2mv.arrayList,G2mv.invarrayList,G2mv.indParmList,
98                    G2mv.invarrayList) = cPickle.load(file)
99                self.parmDict = cPickle.load(file)
100                self.varylist = cPickle.load(file)
101                self.calcControls = cPickle.load(file)
102                self.pawleyLookup = cPickle.load(file)
103                self.use = [False for i in range(len(self.varylist))]
104                self.delt = [max(abs(self.parmDict[name])*0.001,1e-6) for name in self.varylist]
105                file.close()
106                self.UpdateControls(event)
107                print G2mv.VarRemapShow(self.varylist)
108        finally:
109            dlg.Destroy()
110           
111    def UpdateControls(self,event):
112       
113        def OnItemCk(event):
114            Obj = event.GetEventObject()
115            item = ObjInd[Obj.GetId()]
116            self.use[item] = Obj.GetValue()
117           
118        def OnDelValue(event):
119            Obj = event.GetEventObject()
120            item = ObjInd[Obj.GetId()]
121            try:
122                value = float(Obj.GetValue())
123            except ValueError:
124                value = self.delt[item]
125            self.delt[item] = value
126            Obj.SetValue('%g'%(value))
127       
128        self.testDerivPanel.DestroyChildren()
129        ObjInd = {}
130        varylist = self.varylist
131        use = self.use
132        delt = self.delt
133        mainSizer = wx.FlexGridSizer(0,8,5,5)
134        for id,[ck,name,d] in enumerate(zip(use,varylist,delt)):
135            useVal = wx.CheckBox(self.testDerivPanel,label=name)
136            useVal.SetValue(ck)
137            ObjInd[useVal.GetId()] = id
138            useVal.Bind(wx.EVT_CHECKBOX, OnItemCk)
139            mainSizer.Add(useVal,0)
140            delVal = wx.TextCtrl(self.testDerivPanel,wx.ID_ANY,'%g'%(d),style=wx.TE_PROCESS_ENTER)
141            ObjInd[delVal.GetId()] = id
142            delVal.Bind(wx.EVT_TEXT_ENTER,OnDelValue)
143            delVal.Bind(wx.EVT_KILL_FOCUS,OnDelValue)
144            mainSizer.Add(delVal,0)
145#        mainSizer.Layout()
146        self.testDerivPanel.SetSizer(mainSizer)   
147        Size = mainSizer.Fit(self.testDerivPanel)
148        Size[0] = 800
149        Size[1] = max(Size[1],290) + 35
150        self.testDerivPanel.SetScrollbars(10,10,Size[0]/10-4,Size[1]/10-1)
151        self.testDerivPanel.SetSize(Size)
152
153    def OnMakePlots(self,event):
154       
155        def test1():
156            fplot = self.plotNB.add('function test').gca()
157            M = G2stMth.errRefine(self.values,self.HistoPhases,
158                self.parmDict,self.varylist,self.calcControls,
159                self.pawleyLookup,None)
160            fplot.plot(M,'r',label='M')
161            fplot.legend(loc='best')
162           
163        def test2(name,delt):
164            hplot = self.plotNB.add('derivatives test for '+name).gca()
165            dMdV = G2stMth.dervRefine(self.values,self.HistoPhases,self.parmDict,
166                self.varylist,self.calcControls,self.pawleyLookup,None)
167            hplot.plot(dMdV[self.varylist.index(name)],'b',label='analytic deriv')
168            if name in self.varylist:
169                print 'parameter:',name,self.values[self.varylist.index(name)],delt
170                self.values[self.varylist.index(name)] -= delt
171                M0 = G2stMth.errRefine(self.values,self.HistoPhases,self.parmDict,
172                    self.varylist,self.calcControls,self.pawleyLookup,None)
173                self.values[self.varylist.index(name)] += 2.*delt
174                M1 = G2stMth.errRefine(self.values,self.HistoPhases,self.parmDict,
175                    self.varylist,self.calcControls,self.pawleyLookup,None)
176                self.values[self.varylist.index(name)] -= delt   
177                Mn = (M1-M0)/(2.*delt)
178                hplot.plot(Mn,'r',label='numeric deriv')
179                hplot.plot(dMdV[self.varylist.index(name)]-Mn,'g',label='diff')
180            hplot.legend(loc='best')
181           
182        while self.plotNB.nb.GetPageCount():
183            self.plotNB.nb.DeletePage(0)
184        test1()
185        for use,name,delt in zip(self.use,self.varylist,self.delt):
186            if use:
187                test2(name,delt)
188       
189        self.plotNB.Show()
190       
191class testDerivmain(wx.App):
192    def OnInit(self):
193        self.main = testDeriv(None)
194        self.main.Show()
195        self.SetTopWindow(self.main)
196        return True
197
198def main():
199    'Starts main application to compute and plot derivatives'
200    application = testDerivmain(0)
201    application.MainLoop()
202   
203if __name__ == '__main__':
204    main()
Note: See TracBrowser for help on using the repository browser.