source: trunk/GSASIIfpaGUI.py @ 4715

Last change on this file since 4715 was 4594, checked in by vondreele, 5 years ago

remove WACV,wx.ALIGN_CENTER_VERTICAL from all Add for wx.VERTICAL wx.BoxSizers?
fix use of inspect.getargspec - now inspect.getfullargspec for python 3.x in G2dataGUI
Add a File.close in G2img_CBF reader

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