source: trunk/GSASIIfpaGUI.py @ 3580

Last change on this file since 3580 was 3580, checked in by toby, 5 years ago

add tube tails to FPA; add grace export to Publ. plot
  • Property svn:eol-style set to native
File size: 29.8 KB
Line 
1# -*- coding: utf-8 -*-
2########### SVN repository information ###################
3# $Date: $
4# $Author: $
5# $Revision: $
6# $URL: $
7# $Id: $
8########### SVN repository information ###################
9'''
10*GSASIIfpaGUI: Fundamental Parameters Routines*
11===============================================
12
13This module contains routines for getting Fundamental Parameters
14Approach (FPA) input, setting up for running the NIST XRD Fundamental
15Parameters Code, plotting the convolutors and computing a set of peaks
16generated by that code.
17
18'''
19from __future__ import division, print_function
20import wx
21import os.path
22import numpy as np
23
24import NIST_profile as FP
25
26import GSASIIpath
27import GSASIIctrlGUI as G2G
28import GSASIIdataGUI as G2gd
29import GSASIIplot as G2plt
30import GSASIImath as G2mth
31import GSASIIpwd as G2pwd
32
33simParms = {}
34'''Parameters to set range for pattern simulation
35'''
36
37parmDict = {'numWave':2}
38'''Parameter dict used for reading Topas-style values. These are
39converted to SI units and placed into :data:`NISTparms`
40'''
41
42NISTparms = {}
43'''Parameters in a nested dict, with an entry for each concolutor. Entries in
44those dicts have values in SI units (of course). NISTparms can be
45can be input directly or can be from created from :data:`parmDict`
46by :func:`XferFPAsettings`
47'''
48
49BraggBrentanoParms = [
50    ('divergence', 0.5, 'Bragg-Brentano divergence angle (degrees)'),
51    ('soller_angle', 2.0, 'Soller slit axial divergence (degrees)'),
52    ('Rs', 220, 'Diffractometer radius (mm)'),
53    ('filament_length', 12., 'X-ray tube line focus length (mm)'),
54    ('sample_length', 12., 'Illuminated sample length in axial direction (mm)'),
55    ('receiving_slit_length', 12., 'Length of receiving slit in axial direction (mm)'),
56    ('LAC_cm', 0.,'Linear absorption coef. adjusted for packing density (cm-1)'),
57    ('sample_thickness', 1., 'Depth of sample (mm)'),
58    ('convolution_steps', 8, 'Number of Fourier-space bins per two-theta step'),
59    ('tube-tails_width', 0.04,'Tube filament width, in projection at takeoff angle (mm)'),
60    ('tube-tails_L-tail', -1.,'Left-side tube tails width, in projection (mm)'),   
61    ('tube-tails_R-tail', 1.,'Right-side tube tails width, in projection (mm)'),
62    ('tube-tails_rel-I', 0.001,'Tube tails fractional intensity (no units)'),
63    ]
64'''FPA dict entries used in :func:`MakeTopasFPASizer`. Tuple contains
65a dict key, a default value and a description. These are the parameters
66needed for all Bragg Brentano instruments
67'''
68
69BBPointDetector = [
70    ('receiving_slit_width', 0.2, 'Width of receiving slit (mm)'),]
71'''Additional FPA dict entries used in :func:`MakeTopasFPASizer`
72needed for Bragg Brentano instruments with point detectors.
73'''
74
75BBPSDDetector = [
76    ('lpsd_th2_angular_range', 3.0, 'Angular range observed by PSD (degrees 2Theta)'),
77    ('lpsd_equitorial_divergence', 0.1, 'Equatorial divergence of the primary beam (degrees)'),]
78'''Additional FPA dict entries used in :func:`MakeTopasFPASizer`
79needed for Bragg Brentano instruments with linear (1-D) PSD detectors.
80'''
81
82Citation = '''MH Mendenhall, K Mullen && JP Cline. (2015) J. Res. of NIST 120, 223-251. doi:10.6028/jres.120.014.
83'''
84   
85def SetCu2Wave():
86    '''Set the parameters to the two-line Cu K alpha 1+2 spectrum
87    '''
88    parmDict['wave'] = {i:v for i,v in enumerate((1.540596,1.544493))}
89    parmDict['int'] = {i:v for i,v in enumerate((0.653817, 0.346183))}
90    parmDict['lwidth'] = {i:v for i,v in enumerate((0.501844,0.626579))}
91SetCu2Wave() # use these as default
92
93def MakeTopasFPASizer(G2frame,FPdlg,mode,SetButtonStatus):
94    '''Create a GUI with parameters for the NIST XRD Fundamental Parameters Code.
95    Parameter input is modeled after Topas input parameters.
96
97    :param wx.Window FPdlg: Frame or Dialog where GUI will appear
98    :param str mode: either 'BBpoint' or 'BBPSD' for Bragg-Brentano point detector or
99      (linear) position sensitive detector
100    :param dict parmDict: dict to place parameters. If empty, default values from
101      globals BraggBrentanoParms, BBPointDetector & BBPSDDetector will be placed in
102      the array.
103    :returns: a sizer with the GUI controls
104 
105    '''
106    def _onOK(event):
107        XferFPAsettings(parmDict)
108        SetButtonStatus(done=True) # done=True triggers the simulation
109        FPdlg.Destroy()
110    def _onClose(event):
111        SetButtonStatus()
112        FPdlg.Destroy()
113    def _onAddWave(event):
114        parmDict['numWave'] += 1 
115        wx.CallAfter(MakeTopasFPASizer,G2frame,FPdlg,mode,SetButtonStatus)
116    def _onRemWave(event):
117        parmDict['numWave'] -= 1 
118        wx.CallAfter(MakeTopasFPASizer,G2frame,FPdlg,mode,SetButtonStatus)
119    def _onSetCu5Wave(event):
120        parmDict['wave'] = {i:v for i,v in enumerate((1.534753,1.540596,1.541058,1.54441,1.544721))}
121        parmDict['int'] = {i:v for i,v in enumerate((0.0159, 0.5791, 0.0762, 0.2417, 0.0871))}
122        parmDict['lwidth'] = {i:v for i,v in enumerate((3.6854, 0.437, 0.6, 0.52, 0.62))}
123        parmDict['numWave'] = 5
124        wx.CallAfter(MakeTopasFPASizer,G2frame,FPdlg,mode,SetButtonStatus)
125    def _onSetCu2Wave(event):
126        SetCu2Wave()
127        parmDict['numWave'] = 2
128        wx.CallAfter(MakeTopasFPASizer,G2frame,FPdlg,mode,SetButtonStatus)
129    def _onSetPoint(event):
130        wx.CallAfter(MakeTopasFPASizer,G2frame,FPdlg,'BBpoint',SetButtonStatus)
131    def _onSetPSD(event):
132        wx.CallAfter(MakeTopasFPASizer,G2frame,FPdlg,'BBPSD',SetButtonStatus)
133    def PlotTopasFPA(event):
134        XferFPAsettings(parmDict)
135        ttArr = np.arange(max(0.5,
136                            simParms['plotpos']-simParms['calcwid']),
137                            simParms['plotpos']+simParms['calcwid'],
138                            simParms['step'])
139        intArr = np.zeros_like(ttArr)
140        NISTpk = setupFPAcalc()
141        try:
142            center_bin_idx,peakObj = doFPAcalc(
143                NISTpk,ttArr,simParms['plotpos'],simParms['calcwid'],
144                simParms['step'])
145        except Exception as err:
146            msg = "Error computing convolution, revise input"
147            print(msg)
148            print(err)
149            return
150        G2plt.PlotFPAconvolutors(G2frame,NISTpk)
151        pkPts = len(peakObj.peak)
152        pkMax = peakObj.peak.max()
153        startInd = center_bin_idx-(pkPts//2) #this should be the aligned start of the new data
154        # scale peak so max I=10,000 and add into intensity array
155        if startInd < 0:
156            intArr[:startInd+pkPts] += 10000 * peakObj.peak[-startInd:]/pkMax
157        elif startInd > len(intArr):
158            return
159        elif startInd+pkPts >= len(intArr):
160            offset = pkPts - len( intArr[startInd:] )
161            intArr[startInd:startInd+pkPts-offset] += 10000 * peakObj.peak[:-offset]/pkMax
162        else:
163            intArr[startInd:startInd+pkPts] += 10000 * peakObj.peak/pkMax
164        G2plt.PlotXY(G2frame, [(ttArr, intArr)],
165                     labelX=r'$2\theta, deg$',
166                     labelY=r'Intensity (arbitrary)',
167                     Title='FPA peak', newPlot=True, lines=True)
168
169    if FPdlg.GetSizer(): FPdlg.GetSizer().Clear(True)
170    numWave = parmDict['numWave']
171    if mode == 'BBpoint':
172        itemList = BraggBrentanoParms+BBPointDetector
173    elif mode == 'BBPSD':
174        itemList = BraggBrentanoParms+BBPSDDetector
175    else:
176        raise Exception('Unknown mode in MakeTopasFPASizer: '+mode)
177   
178    MainSizer = wx.BoxSizer(wx.VERTICAL)
179    MainSizer.Add((-1,5))
180    waveSizer = wx.FlexGridSizer(cols=numWave+1,hgap=3,vgap=5)
181    for lbl,prm,defVal in zip(
182            (u'Wavelength (\u212b)','Rel. Intensity',u'Lorentz Width\n(\u212b/1000)'),
183            ('wave','int','lwidth'),
184            (0.0,   1.0,   0.1),
185            ):
186        text = wx.StaticText(FPdlg,wx.ID_ANY,lbl,style=wx.ALIGN_CENTER)
187        text.SetBackgroundColour(wx.WHITE)
188        waveSizer.Add(text,0,wx.EXPAND)
189        if prm not in parmDict: parmDict[prm] = {}
190        for i in range(numWave):
191            if i not in parmDict[prm]: parmDict[prm][i] = defVal
192            ctrl = G2G.ValidatedTxtCtrl(FPdlg,parmDict[prm],i,size=(90,-1))
193            waveSizer.Add(ctrl,1,wx.ALIGN_CENTER_VERTICAL,1)
194    MainSizer.Add(waveSizer)
195    MainSizer.Add((-1,5))
196    btnsizer = wx.BoxSizer(wx.HORIZONTAL)
197    btn = wx.Button(FPdlg, wx.ID_ANY,'Add col')
198    btnsizer.Add(btn)
199    btn.Bind(wx.EVT_BUTTON,_onAddWave)
200    btn = wx.Button(FPdlg, wx.ID_ANY,'Remove col')
201    btnsizer.Add(btn)
202    btn.Bind(wx.EVT_BUTTON,_onRemWave)
203    btn = wx.Button(FPdlg, wx.ID_ANY,'CuKa1+2')
204    btnsizer.Add(btn)
205    btn.Bind(wx.EVT_BUTTON,_onSetCu2Wave)
206    btn = wx.Button(FPdlg, wx.ID_ANY,'CuKa-5wave')
207    btnsizer.Add(btn)
208    btn.Bind(wx.EVT_BUTTON,_onSetCu5Wave)
209    MainSizer.Add(btnsizer, 0, wx.ALIGN_CENTER, 0)
210    MainSizer.Add((-1,5))
211    btnsizer = wx.BoxSizer(wx.HORIZONTAL)
212    btn = wx.Button(FPdlg, wx.ID_ANY,'Point Dect.')
213    btn.Enable(not mode == 'BBpoint')
214    btnsizer.Add(btn)
215    btn.Bind(wx.EVT_BUTTON,_onSetPoint)
216    btn = wx.Button(FPdlg, wx.ID_ANY,'PSD')
217    btn.Enable(not mode == 'BBPSD')
218    btnsizer.Add(btn)
219    btn.Bind(wx.EVT_BUTTON,_onSetPSD)
220    MainSizer.Add(btnsizer, 0, wx.ALIGN_CENTER, 0)
221    MainSizer.Add((-1,5))
222   
223    prmSizer = wx.FlexGridSizer(cols=3,hgap=3,vgap=5)
224    text = wx.StaticText(FPdlg,wx.ID_ANY,'label',style=wx.ALIGN_CENTER)
225    text.SetBackgroundColour(wx.WHITE)
226    prmSizer.Add(text,0,wx.EXPAND)
227    text = wx.StaticText(FPdlg,wx.ID_ANY,'value',style=wx.ALIGN_CENTER)
228    text.SetBackgroundColour(wx.WHITE)
229    prmSizer.Add(text,0,wx.EXPAND)
230    text = wx.StaticText(FPdlg,wx.ID_ANY,'explanation',style=wx.ALIGN_CENTER)
231    text.SetBackgroundColour(wx.WHITE)
232    prmSizer.Add(text,0,wx.EXPAND)
233    for lbl,defVal,text in itemList:
234        prmSizer.Add(wx.StaticText(FPdlg,wx.ID_ANY,lbl),1,wx.ALIGN_RIGHT|wx.ALIGN_CENTER_VERTICAL,1)
235        if lbl not in parmDict: parmDict[lbl] = defVal
236        ctrl = G2G.ValidatedTxtCtrl(FPdlg,parmDict,lbl,size=(70,-1))
237        prmSizer.Add(ctrl,1,wx.ALL|wx.ALIGN_CENTER_VERTICAL,1)
238        txt = wx.StaticText(FPdlg,wx.ID_ANY,text,size=(400,-1))
239        txt.Wrap(380)
240        prmSizer.Add(txt)
241    MainSizer.Add(prmSizer)
242    MainSizer.Add((-1,4),1,wx.EXPAND,1)
243    btnsizer = wx.BoxSizer(wx.HORIZONTAL)
244    btn = wx.Button(FPdlg, wx.ID_ANY, 'Plot peak')
245    btnsizer.Add(btn)
246    btn.Bind(wx.EVT_BUTTON,PlotTopasFPA)
247    btnsizer.Add(wx.StaticText(FPdlg,wx.ID_ANY,' at '))
248    if 'plotpos' not in simParms: simParms['plotpos'] =  simParms['minTT']
249    ctrl = G2G.ValidatedTxtCtrl(FPdlg,simParms,'plotpos',size=(70,-1))
250    btnsizer.Add(ctrl)
251    btnsizer.Add(wx.StaticText(FPdlg,wx.ID_ANY,' deg.'))   
252    MainSizer.Add(btnsizer, 0, wx.ALIGN_CENTER, 0)
253    MainSizer.Add((-1,4),1,wx.EXPAND,1)
254    btnsizer = wx.BoxSizer(wx.HORIZONTAL)
255    OKbtn = wx.Button(FPdlg, wx.ID_OK)
256    OKbtn.SetDefault()
257    btnsizer.Add(OKbtn)
258    Cbtn = wx.Button(FPdlg, wx.ID_CLOSE,"Cancel") 
259    btnsizer.Add(Cbtn)
260    MainSizer.Add(btnsizer, 0, wx.ALIGN_CENTER, 0)
261    MainSizer.Add((-1,4),1,wx.EXPAND,1)
262    # bindings for close of window
263    OKbtn.Bind(wx.EVT_BUTTON,_onOK)
264    Cbtn.Bind(wx.EVT_BUTTON,_onClose)
265    FPdlg.SetSizer(MainSizer)
266    MainSizer.Layout()
267    MainSizer.Fit(FPdlg)
268    FPdlg.SetMinSize(FPdlg.GetSize())
269    FPdlg.SendSizeEvent()
270
271def XferFPAsettings(InpParms):
272    '''convert Topas-type parameters to SI units for NIST and place in a dict sorted
273    according to use in each convoluter
274
275    :param dict InpParms: a dict with Topas-like parameters, as set in
276      :func:`MakeTopasFPASizer`
277    :returns: a nested dict with global parameters and those for each convolution
278    '''
279    wavenums = range(InpParms['numWave'])
280    source_wavelengths_m = 1.e-10 * np.array([InpParms['wave'][i] for i in wavenums])
281    la = [InpParms['int'][i] for i in wavenums]
282    source_intensities = np.array(la)/max(la)
283    source_lor_widths_m = 1.e-10 * 1.e-3 * np.array([InpParms['lwidth'][i] for i in wavenums])
284    source_gauss_widths_m = 1.e-10 * 1.e-3 * np.array([0.001 for i in wavenums])
285   
286    NISTparms["emission"] = {'emiss_wavelengths' : source_wavelengths_m,
287                'emiss_intensities' : source_intensities,
288                'emiss_gauss_widths' : source_gauss_widths_m,
289                'emiss_lor_widths' : source_lor_widths_m,
290                'crystallite_size_gauss' : 1.e-9 * InpParms.get('Size_G',1e6),
291                'crystallite_size_lor' : 1.e-9 * InpParms.get('Size_L',1e6)}
292   
293    if InpParms['filament_length'] == InpParms['receiving_slit_length']: # workaround:
294        InpParms['receiving_slit_length'] *= 1.00001 # avoid bug when slit lengths are identical
295    NISTparms["axial"]  = {
296            'axDiv':"full", 'slit_length_source' : 1e-3*InpParms['filament_length'],
297            'slit_length_target' : 1e-3*InpParms['receiving_slit_length'],
298            'length_sample' : 1e-3 * InpParms['sample_length'], 
299            'n_integral_points' : 10,
300            'angI_deg' : InpParms['soller_angle'],
301            'angD_deg': InpParms['soller_angle']
302            }
303    if InpParms.get('LAC_cm',0) > 0:
304        NISTparms["absorption"] = {
305            'absorption_coefficient': InpParms['LAC_cm']*100, #like LaB6, in m^(-1)
306            'sample_thickness': 1e-3 * InpParms['sample_thickness'],
307            }
308    elif "absorption" in NISTparms:
309        del NISTparms["absorption"]
310
311    if InpParms.get('lpsd_equitorial_divergence',0) > 0 and InpParms.get(
312            'lpsd_th2_angular_range',0) > 0:
313        PSDdetector_length_mm=np.arcsin(np.pi*InpParms['lpsd_th2_angular_range']/180.
314                                            )*InpParms['Rs'] # mm
315        NISTparms["si_psd"] = {
316            'equatorial_divergence_deg': InpParms['lpsd_equitorial_divergence'],
317            'si_psd_window_bounds': (0.,PSDdetector_length_mm/1000.)
318            }
319    elif "si_psd" in NISTparms:
320        del NISTparms["si_psd"]
321       
322    if InpParms.get('Specimen_Displacement'):
323        NISTparms["displacement"] = {'specimen_displacement': 1e-3 * InpParms['Specimen_Displacement']}
324    elif "displacement" in NISTparms:
325        del NISTparms["displacement"]
326
327    if InpParms.get('receiving_slit_width'):
328        NISTparms["receiver_slit"] = {'slit_width':1e-3*InpParms['receiving_slit_width']}
329    elif "receiver_slit" in NISTparms:
330        del NISTparms["receiver_slit"]
331
332    if InpParms.get('tube-tails_width', 0) > 0 and InpParms.get(
333            'tube-tails_rel-I',0) > 0:
334        NISTparms["tube_tails"] = {
335            'main_width' : 1e-3 * InpParms.get('tube-tails_width', 0.),
336            'tail_left' : -1e-3 * InpParms.get('tube-tails_L-tail',0.),
337            'tail_right' : 1e-3 * InpParms.get('tube-tails_R-tail',0.),
338            'tail_intens' : InpParms.get('tube-tails_rel-I',0.),}
339    elif "tube_tails" in NISTparms:
340        del NISTparms["tube_tails"]
341
342    # set Global parameters
343    max_wavelength = source_wavelengths_m[np.argmax(source_intensities)]
344    NISTparms[""] = {
345        'equatorial_divergence_deg' : InpParms['divergence'],
346        'dominant_wavelength' : max_wavelength,
347        'diffractometer_radius' : 1e-3* InpParms['Rs'],
348        'oversampling' : InpParms['convolution_steps'],
349        }
350def setupFPAcalc():
351    '''Create a peak profile object using the NIST XRD Fundamental
352    Parameters Code.
353   
354    :returns: a profile object that can provide information on
355      each convolution or compute the composite peak shape.
356    '''
357    p=FP.FP_profile(anglemode="twotheta",
358                    output_gaussian_smoother_bins_sigma=1.0,
359                    oversampling=NISTparms.get('oversampling',10))
360    p.debug_cache=False
361    #set parameters for each convolver
362    for key in NISTparms:
363        if key:
364            p.set_parameters(convolver=key,**NISTparms[key])
365        else:
366            p.set_parameters(**NISTparms[key])
367    return p
368       
369def doFPAcalc(NISTpk,ttArr,twotheta,calcwid,step):
370    '''Compute a single peak using a NIST profile object
371
372    :param object NISTpk: a peak profile computational object from the
373      NIST XRD Fundamental Parameters Code, typically established from
374      a call to :func:`SetupFPAcalc`
375    :param np.Array ttArr: an evenly-spaced grid of two-theta points (degrees)
376    :param float twotheta: nominal center of peak (degrees)
377    :param float calcwid: width to perform convolution (degrees)
378    :param float step: step size
379    '''
380    # find closest point to twotheta (may be outside limits of the array)
381    center_bin_idx=min(ttArr.searchsorted(twotheta),len(ttArr)-1)
382    NISTpk.set_optimized_window(twotheta_exact_bin_spacing_deg=step,
383                twotheta_window_center_deg=ttArr[center_bin_idx],
384                twotheta_approx_window_fullwidth_deg=calcwid,
385                )
386    NISTpk.set_parameters(twotheta0_deg=twotheta)
387    return center_bin_idx,NISTpk.compute_line_profile()
388
389def MakeSimSizer(G2frame, dlg):
390    '''Create a GUI to get simulation with parameters for Fundamental
391    Parameters fitting.
392
393    :param wx.Window dlg: Frame or Dialog where GUI will appear
394
395    :returns: a sizer with the GUI controls
396 
397    '''
398    def _onOK(event):
399        msg = ''
400        if simParms['minTT']-simParms['calcwid']/1.5 < 0.1:
401            msg += 'First peak minus half the calc width is too low'
402        if simParms['maxTT']+simParms['calcwid']/1.5 > 175:
403            if msg: msg += '\n'
404            msg += 'Last peak plus half the calc width is too high'
405        if simParms['npeaks'] < 8:
406            if msg: msg += '\n'
407            msg += 'At least 8 peaks are needed'
408        if msg:
409            G2G.G2MessageBox(dlg,msg,'Bad input, try again')
410            return
411        # compute "obs" pattern
412        ttArr = np.arange(max(0.5,
413                            simParms['minTT']-simParms['calcwid']/1.5),
414                            simParms['maxTT']+simParms['calcwid']/1.5,
415                            simParms['step'])
416        intArr = np.zeros_like(ttArr)
417        peaklist = np.linspace(simParms['minTT'],simParms['maxTT'],
418                               simParms['npeaks'],endpoint=True)
419        peakSpacing = (peaklist[-1]-peaklist[0])/(len(peaklist)-1)
420        NISTpk = setupFPAcalc()
421        minPtsHM = len(intArr)  # initialize points above half-max
422        maxPtsHM = 0
423        for num,twoth_peak in enumerate(peaklist):
424            try:
425                center_bin_idx,peakObj = doFPAcalc(
426                    NISTpk,ttArr,twoth_peak,simParms['calcwid'],
427                    simParms['step'])
428            except:
429                if msg: msg += '\n'
430                msg = "Error computing convolution, revise input"
431                continue
432            if num == 0: G2plt.PlotFPAconvolutors(G2frame,NISTpk)
433            pkMax = peakObj.peak.max()
434            pkPts = len(peakObj.peak)
435            minPtsHM = min(minPtsHM,sum(peakObj.peak >= 0.5*pkMax)) # points above half-max
436            maxPtsHM = max(maxPtsHM,sum(peakObj.peak >= 0.5*pkMax)) # points above half-max
437            startInd = center_bin_idx-(pkPts//2) #this should be the aligned start of the new data
438            # scale peak so max I=10,000 and add into intensity array
439            if startInd < 0:
440                intArr[:startInd+pkPts] += 10000 * peakObj.peak[-startInd:]/pkMax
441            elif startInd > len(intArr):
442                break
443            elif startInd+pkPts >= len(intArr):
444                offset = pkPts - len( intArr[startInd:] )
445                intArr[startInd:startInd+pkPts-offset] += 10000 * peakObj.peak[:-offset]/pkMax
446            else:
447                intArr[startInd:startInd+pkPts] += 10000 * peakObj.peak/pkMax
448        # check if peaks are too closely spaced
449        if maxPtsHM*simParms['step'] > peakSpacing/4:
450            if msg: msg += '\n'
451            msg += 'Maximum FWHM ({}) is too large compared to the peak spacing ({}). Decrease number of peaks or increase data range.'.format(
452                maxPtsHM*simParms['step'], peakSpacing)
453        # check if too few points across Hmax
454        if minPtsHM < 10:
455            if msg: msg += '\n'
456            msg += 'There are only {} points above the half-max. 10 are needed. Dropping step size.'.format(minPtsHM)
457            simParms['step'] *= 0.5
458        if msg:
459            G2G.G2MessageBox(dlg,msg,'Bad input, try again')
460            wx.CallAfter(MakeSimSizer,G2frame, dlg)
461            return
462        # pattern has been computed successfully
463        dlg.Destroy()
464        wx.CallAfter(FitFPApeaks,ttArr, intArr, peaklist, maxPtsHM) # do peakfit outside event callback
465
466    def FitFPApeaks(ttArr, intArr, peaklist, maxPtsHM):
467        '''Perform a peak fit to the FP simulated pattern
468        '''
469        plswait = wx.Dialog(G2frame,style=wx.DEFAULT_DIALOG_STYLE | wx.RESIZE_BORDER)
470        vbox = wx.BoxSizer(wx.VERTICAL)
471        vbox.Add((1,1),1,wx.ALL|wx.EXPAND,1)
472        txt = wx.StaticText(plswait,wx.ID_ANY,
473                                'Fitting peaks...\nPlease wait...',
474                                style=wx.ALIGN_CENTER)
475        vbox.Add(txt,0,wx.ALL|wx.EXPAND)
476        vbox.Add((1,1),1,wx.ALL|wx.EXPAND,1)
477        plswait.SetSizer(vbox)
478        plswait.Layout()
479        plswait.CenterOnParent()
480        plswait.Show() # post "please wait"
481        wx.BeginBusyCursor()
482        # pick out one or two most intense wavelengths
483        ints = list(NISTparms['emission']['emiss_intensities'])
484        Lam1 = NISTparms['emission']['emiss_wavelengths'][np.argmax(ints)]*1e10
485        if len(ints) > 1: 
486            ints[np.argmax(ints)] = -1
487            Lam2 = NISTparms['emission']['emiss_wavelengths'][np.argmax(ints)]*1e10
488        else:
489            Lam2 = None
490        histId = G2frame.AddSimulatedPowder(ttArr,intArr,
491                                       'NIST Fundamental Parameters simulation',
492                                       Lam1,Lam2)
493        controls = G2frame.GPXtree.GetItemPyData(
494            G2gd.GetGPXtreeItemId(G2frame,G2frame.root,'Controls'))
495        controldat = controls.get('data',
496                            {'deriv type':'analytic','min dM/M':0.001,})  #fil
497        Parms,Parms2 = G2frame.GPXtree.GetItemPyData(
498            G2gd.GetGPXtreeItemId(G2frame,histId,'Instrument Parameters'))
499        peakData = G2frame.GPXtree.GetItemPyData(
500            G2gd.GetGPXtreeItemId(G2frame,histId,'Peak List'))
501        # set background to 0 with one term = 0; disable refinement
502        bkg1,bkg2 = bkg = G2frame.GPXtree.GetItemPyData(
503            G2gd.GetGPXtreeItemId(G2frame,histId,'Background'))
504        bkg1[1]=False
505        bkg1[2]=0
506        bkg1[3]=0.0
507        limits = G2frame.GPXtree.GetItemPyData(
508            G2gd.GetGPXtreeItemId(G2frame,histId,'Limits'))
509        # approximate asym correction
510        try:
511            Parms['SH/L'][1] = 0.25 * (
512                NISTparms['axial']['length_sample']+
513                NISTparms['axial']['slit_length_source']
514                        ) / NISTparms['']['diffractometer_radius']
515        except:
516            pass
517           
518        for pos in peaklist:
519            i = ttArr.searchsorted(pos)
520            area = sum(intArr[max(0,i-maxPtsHM):min(len(intArr),i+maxPtsHM)])
521            peakData['peaks'].append(G2mth.setPeakparms(Parms,Parms2,pos,area))
522        histData = G2frame.GPXtree.GetItemPyData(histId)
523        # refine peak positions only
524        bxye = np.zeros(len(histData[1][1]))
525        peakData['sigDict'] = G2pwd.DoPeakFit('LSQ',peakData['peaks'],
526                                            bkg,limits[1],
527                                            Parms,Parms2,histData[1],bxye,[],
528                                           False,controldat,None)[0]
529        # refine peak areas as well
530        for pk in peakData['peaks']:
531            pk[1] = True
532        peakData['sigDict'] = G2pwd.DoPeakFit('LSQ',peakData['peaks'],
533                                            bkg,limits[1],
534                                            Parms,Parms2,histData[1],bxye,[],
535                                           False,controldat)[0]
536        # refine profile function
537        for p in ('U', 'V', 'W', 'X', 'Y'):
538            Parms[p][2] = True
539        peakData['sigDict'] = G2pwd.DoPeakFit('LSQ',peakData['peaks'],
540                                            bkg,limits[1],
541                                            Parms,Parms2,histData[1],bxye,[],
542                                           False,controldat)[0]
543        # add in asymmetry
544        Parms['SH/L'][2] = True
545        peakData['sigDict'] = G2pwd.DoPeakFit('LSQ',peakData['peaks'],
546                                            bkg,limits[1],
547                                            Parms,Parms2,histData[1],bxye,[],
548                                           False,controldat)[0]
549        # reset "initial" profile
550        for p in Parms:
551            if len(Parms[p]) == 3:
552                Parms[p][0] = Parms[p][1]
553                Parms[p][2] = False
554        wx.EndBusyCursor()
555        plswait.Destroy() # remove "please wait"
556        # save Iparms
557        pth = G2G.GetExportPath(G2frame)
558        fldlg = wx.FileDialog(G2frame, 'Set name to save GSAS-II instrument parameters file', pth, '', 
559            'instrument parameter files (*.instprm)|*.instprm',wx.FD_SAVE|wx.FD_OVERWRITE_PROMPT)
560        try:
561            if fldlg.ShowModal() == wx.ID_OK:
562                filename = fldlg.GetPath()
563                # make sure extension is .instprm
564                filename = os.path.splitext(filename)[0]+'.instprm'
565                File = open(filename,'w')
566                File.write("#GSAS-II instrument parameter file; do not add/delete items!\n")
567                for item in Parms:
568                    File.write(item+':'+str(Parms[item][1])+'\n')
569                File.close()
570                print ('Instrument parameters saved to: '+filename)
571        finally:
572            fldlg.Destroy()
573        #GSASIIpath.IPyBreak()
574       
575    def _onClose(event):
576        dlg.Destroy()
577    def SetButtonStatus(done=False):
578        OKbtn.Enable(bool(NISTparms))
579        saveBtn.Enable(bool(NISTparms))
580        if done: _onOK(None)
581    def _onSetFPA(event):
582        # Create a non-modal dialog for Topas-style FP input.
583        FPdlg = wx.Dialog(dlg,wx.ID_ANY,'FPA parameters',
584                style=wx.DEFAULT_DIALOG_STYLE|wx.RESIZE_BORDER)
585        MakeTopasFPASizer(G2frame,FPdlg,'BBpoint',SetButtonStatus)
586        FPdlg.CenterOnParent()
587        FPdlg.Raise()
588        FPdlg.Show() 
589    def _onSaveFPA(event):
590        filename = G2G.askSaveFile(G2frame,'','.NISTfpa',
591                                       'dict of NIST FPA values',dlg)
592        if not filename: return
593        fp = open(filename,'w')
594        fp.write('# parameters to be used in the NIST XRD Fundamental Parameters program\n')
595        fp.write('{\n')
596        for key in sorted(NISTparms):
597            fp.write("   '"+key+"' : "+str(NISTparms[key])+",")
598            if not key: fp.write('  # global parameters')
599            fp.write('\n')
600        fp.write('}\n')
601        fp.close()
602    def _onReadFPA(event):
603        filename = G2G.GetImportFile(G2frame,
604                message='Read file with dict of values for NIST Fundamental Parameters',
605                parent=dlg,
606                wildcard='dict of NIST FPA values|*.NISTfpa')
607        if not filename: return
608        if not filename[0]: return
609        try:
610            txt = open(filename[0],'r').read()
611            NISTparms.clear()
612            array = np.array
613            d = eval(txt)
614            NISTparms.update(d)
615        except Exception as err:
616            G2G.G2MessageBox(dlg,
617                    u'Error reading file {}:{}\n'.format(filename,err),
618                    'Bad dict input')
619        #GSASIIpath.IPyBreak()
620        SetButtonStatus()
621
622    if dlg.GetSizer(): dlg.GetSizer().Clear(True)
623    MainSizer = wx.BoxSizer(wx.VERTICAL)
624    MainSizer.Add(wx.StaticText(dlg,wx.ID_ANY,
625            'Fit Profile Parameters to Peaks from Fundamental Parameters',
626            style=wx.ALIGN_CENTER),0,wx.EXPAND)
627    MainSizer.Add((-1,5))
628    prmSizer = wx.FlexGridSizer(cols=2,hgap=3,vgap=5)
629    text = wx.StaticText(dlg,wx.ID_ANY,'value',style=wx.ALIGN_CENTER)
630    text.SetBackgroundColour(wx.WHITE)
631    prmSizer.Add(text,0,wx.EXPAND)
632    text = wx.StaticText(dlg,wx.ID_ANY,'explanation',style=wx.ALIGN_CENTER)
633    text.SetBackgroundColour(wx.WHITE)
634    prmSizer.Add(text,0,wx.EXPAND)
635    for key,defVal,text in (
636            ('minTT',3.,'Location of first peak in 2theta (deg)'),
637            ('maxTT',123.,'Location of last peak in 2theta (deg)'),
638            ('step',0.01,'Pattern step size (deg 2theta)'),
639            ('npeaks',13.,'Number of peaks'),
640            ('calcwid',2.,'Range to compute each peak (deg 2theta)'),
641            ):
642        if key not in simParms: simParms[key] = defVal
643        ctrl = G2G.ValidatedTxtCtrl(dlg,simParms,key,size=(70,-1))
644        prmSizer.Add(ctrl,1,wx.ALL|wx.ALIGN_CENTER_VERTICAL,1)
645        txt = wx.StaticText(dlg,wx.ID_ANY,text,size=(300,-1))
646        txt.Wrap(280)
647        prmSizer.Add(txt)
648    MainSizer.Add(prmSizer)
649
650    btnsizer = wx.BoxSizer(wx.HORIZONTAL)
651    btn = wx.Button(dlg, wx.ID_ANY,'Input FP vals')
652    btnsizer.Add(btn)
653    btn.Bind(wx.EVT_BUTTON,_onSetFPA)
654    saveBtn = wx.Button(dlg, wx.ID_ANY,'Save FPA dict')
655    btnsizer.Add(saveBtn)
656    saveBtn.Bind(wx.EVT_BUTTON,_onSaveFPA)
657    readBtn = wx.Button(dlg, wx.ID_ANY,'Read FPA dict')
658    btnsizer.Add(readBtn)
659    readBtn.Bind(wx.EVT_BUTTON,_onReadFPA)
660    MainSizer.Add(btnsizer, 0, wx.ALIGN_CENTER, 0)
661    MainSizer.Add((-1,4),1,wx.EXPAND,1)
662    txt = wx.StaticText(dlg,wx.ID_ANY,
663                            'If you use this, please cite: '+Citation,
664                            size=(350,-1))
665    txt.Wrap(340)
666    MainSizer.Add(txt,0,wx.ALIGN_CENTER)
667    btnsizer = wx.BoxSizer(wx.HORIZONTAL)
668    OKbtn = wx.Button(dlg, wx.ID_OK)
669    OKbtn.SetDefault()
670    btnsizer.Add(OKbtn)
671    Cbtn = wx.Button(dlg, wx.ID_CLOSE,"Cancel") 
672    btnsizer.Add(Cbtn)
673    MainSizer.Add(btnsizer, 0, wx.ALIGN_CENTER, 0)
674    MainSizer.Add((-1,4),1,wx.EXPAND,1)
675    # bindings for close of window
676    OKbtn.Bind(wx.EVT_BUTTON,_onOK)
677    Cbtn.Bind(wx.EVT_BUTTON,_onClose)
678    SetButtonStatus()
679    dlg.SetSizer(MainSizer)
680    MainSizer.Layout()
681    MainSizer.Fit(dlg)
682    dlg.SetMinSize(dlg.GetSize())
683    dlg.SendSizeEvent()
684    dlg.Raise()
685   
686def GetFPAInput(G2frame):
687    dlg = wx.Dialog(G2frame,wx.ID_ANY,'FPA input',
688            style=wx.DEFAULT_DIALOG_STYLE|wx.RESIZE_BORDER)
689    MakeSimSizer(G2frame,dlg)
690    dlg.CenterOnParent()
691    dlg.Show()
692    return
693       
Note: See TracBrowser for help on using the repository browser.