source: Tutorials/PythonScript/Scripting.htm @ 4131

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

add fit profile terms tutorial

  • Property svn:mime-type set to text/html
File size: 28.9 KB
Line 
1<!DOCTYPE HTML PUBLIC "-//IETF//DTD HTML//EN">
2<html> <head>
3<title>Scripting GSAS-II</title>
4</head>
5
6<body>
7<h1>Running a GSAS-II Refinement From a Python Script</h1>
8To demonstrate the use of the
9<A href="http://gsas-ii.readthedocs.io/en/latest/GSASIIscripts.html" target="_blank">
10GSASIIscriptable module</A>. This uses a Python script
11to perform a refinement or computation, but without use of the GSAS-II
12graphical user interface. Note that the
13<A href="http://gsas-ii.readthedocs.io/en/latest/GSASIIscripts.html" target="_blank">
14GSASIIscriptable module</A> does not offer access to all of the
15features within GSAS-II, but over time it is expected to grow in
16capabilities. Much of the initial development for this module was done
17by Jackson O'Donnell, as a summer undergraduate visitor under
18supervisor Dr. Maria Chan. Other programming contributions are welcome.
19<P>
20This tutorial has three sections:
21<OL type="A">
22<LI><a href="#multistep"> Create a Multi-step Script</A></LI>
23<LI><a href="#SingleStep">Single-Step Refinement</a></LI>
24<LI><a href="#Simulate">Powder Pattern Simulation</a></LI>
25</OL>
26The first section duplicates the refinement in the
27<A
28href="https://subversion.xray.aps.anl.gov/pyGSAS/Tutorials/CWCombined/Combined%20refinement.htm"
29target="_blank">GSAS-II CW Combined Refinement</A>
30tutorial as a multi-step process. The second section repeats the previous
31refinement, but demonstrates how a complex process can be entered into
32a single Python dict and the refinement executed with a single
33function call. The third section shows how a pattern simulation,
34rather than refinement can be executed from a script.
35
36
37<h2>Prerequisites</h2>
38This exercise assumes that the reader has reasonable familiarity with
39the Python language.
40Also, the documentation on the
41<A href="http://gsas-ii.readthedocs.io/en/latest/GSASIIscripts.html"
42target="_blank">
43Scripting Tools in the GSAS-II GSASIIscriptable module</A> should be
44read from here: <A href="http://gsas-ii.readthedocs.io/en/latest/GSASIIscripts.html"
45target="_blank">http://gsas-ii.readthedocs.io/en/latest/GSASIIscripts.html</A>. It
46is also wise to read the <A href="https://subversion.xray.aps.anl.gov/pyGSAS/Tutorials/CWCombined/Combined%20refinement.htm"
47target="_blank">GSAS-II CW Combined Refinement</A> tutorial that this
48exercise is modeled upon, which explains why each refinement step is
49being used.
50The exercise will require downloading of the same files needed for the
51Combined Refinement tutorial: "PBSO4.XRA", "INST_XRY.PRM", "PBSO4.CWN",
52"inst_d1a.prm" and "PbSO4-Wyckoff.cif",
53which can be downloaded from
54<A href="https://subversion.xray.aps.anl.gov/pyGSAS/Tutorials/PythonScript/data/">
55https://subversion.xray.aps.anl.gov/pyGSAS/Tutorials/PythonScript/data/</A>)
56or will be downloaded with this tutorial if that is requested.
57<P>
58The exercise can be performed by placing all of the Python commands
59into a script,
60but a more pedagogical approach will be to enter the
61commands into a Python interpreter. Use of IPython or Jupyter to run
62Python will make this a more pleasant experience.
63
64<a name="multistep">
65<h2>A. Create a Multi-step Script</h2></a>
66
67Note that the script in this section is supplied for
68<A href="https://subversion.xray.aps.anl.gov/trac/pyGSAS/browser/Tutorials/PythonScript/data/example.py?format=txt">
69download here</A>, but some editing will be needed to match where files
70are found on your computer.
71
72<h4>0: Load the GSASIIscriptable module</H4>
73<blockquote>
74The first step in writing any Python module is to load the modules
75that will be needed. Here we will need the <tt>os</tt> and
76<tt>sys</tt> modules from the Python standard library and the
77GSASIIscriptable from the GSAS-II code. The location where the
78GSAS-II source code is installed must usually be hard-coded into the
79script as is done in the example below. Note that a common location
80for this will be
81<TT>os.path.expanduser("~/g2conda/GSASII/")</TT> or
82<tt>'/Users/toby/software/GSASII'</tt>, etc.
83Thus the script will begin with something like this:
84
85<blockquote><textarea rows="3" cols="60" readonly>
86import os,sys
87sys.path.insert(0,os.path.expanduser("~/g2conda/GSASII/"))
88import GSASIIscriptable as G2sc</textarea></blockquote>
89
90If a <tt>ImportError: No module named GSASIIscriptable</tt> error
91occurs, then the wrong directory location has been supplied for the
92GSAS-II files.
93</blockquote>
94
95<h4>1: Define some other prerequisite code</H4>
96<blockquote>
97To simplify this example, we will define the location where files will
98be written as <tt>workdir</tt> (this directory must exist, but it may
99be empty) and the location where the input files for this exercise
100("PBSO4.XRA", "INST_XRY.PRM", "PBSO4.CWN",
101"inst_d1a.prm" and "PbSO4-Wyckoff.cif") will
102be found as <tt>datadir</tt>. (As discussed previously,
103these files can be
104<A href="https://subversion.xray.aps.anl.gov/pyGSAS/Tutorials/PythonScript/data/">
105downloaded from here. </A>)
106<P>
107We also define here a short function to display the weighted profile R
108factor for every histogram in a project, <tt>HistStats</tt>. This will
109be <A href="#HistStats">discussed later</A>.
110
111<blockquote><textarea rows="10" cols="70" readonly>
112workdir = "/Users/toby/Scratch/PythonScript"
113datadir = "/Users/toby/software/G2/Tutorials/PythonScript/data"
114
115def HistStats(gpx):
116    '''prints profile rfactors for all histograms'''
117    print(u"*** profile Rwp, "+os.path.split(gpx.filename)[1])
118    for hist in gpx.histograms():
119        print("\t{:20s}: {:.2f}".format(hist.name,hist.get_wR()))
120    print("")
121    gpx.save()</textarea></blockquote>
122
123</blockquote>
124
125<h4>2: Create a GSAS-II project</H4>
126<blockquote>
127The first step in creating a GSASIIscriptable script to to create or
128access a GSAS-II project, which is done with a call to
129<TT>GSASIIscriptable.G2Project()</tt>. This can be done in one of two
130ways: a call with <tt>G2sc.G2Project(newgpx=</tt><I>file</I><tt>)</tt> creates a new
131(empty) project, while a call with
132<tt>G2sc.G2Project(gpxfile=</tt><I>file</I><tt>)</tt> opens and
133reads an existing project (.gpx) file. Both return a <tt>G2Project</tt> wrapper
134object that is used to access a number of methods and variables.
135Note that <TT>GSASIIscriptable.G2Project()</tt> can read .gpx files
136written by both previous <TT>GSASIIscriptable.G2Project()</tt> runs or
137the GSAS-II GUI.
138<P>
139In this example, we create a new project by using the
140<tt>newgpx=</tt><I>file</I> argument with this code:
141   
142<blockquote><textarea rows="3" cols="70" readonly>
143# create a project with a default project name
144gpx = G2sc.G2Project(newgpx='PbSO4.gpx')</textarea></blockquote>
145
146Note that the file will not actually be created until an operation
147that saves the project is called.
148</blockquote>
149
150<h4>3: Add Histograms and Phase to the GSAS-II project</H4>
151<blockquote>
152To add two powder diffraction datasets (histograms) and a phase to the
153project using this code:
154
155<blockquote><textarea rows="10" cols="70" readonly>
156# setup step 1: add two histograms to the project
157hist1 = gpx.add_powder_histogram(os.path.join(datadir,"PBSO4.XRA"),
158                          os.path.join(datadir,"INST_XRY.PRM"))
159hist2 = gpx.add_powder_histogram(os.path.join(datadir,"PBSO4.CWN"),
160                          os.path.join(datadir,"inst_d1a.prm"))
161# setup step 2: add a phase and link it to the previous histograms
162phase0 = gpx.add_phase(os.path.join(datadir,"PbSO4-Wyckoff.cif"),
163                      phasename="PbSO4",
164                      histograms=[hist1,hist2])</textarea></blockquote>
165
166 
167We use the <I>project</I><tt>.add_powder_histogram()</tt> method in
168the <TT>GSASIIscriptable</tt> class to read in the powder data.
169The two arguments to
170<tt>.add_powder_histogram()</tt> are the powder dataset and the
171  instrument parameter file. Note that these files are both "GSAS
172  powder data" files and the importer for this format (and many
173  others) allows files with arbitrary extensions to be read.
174All importers that allow for extensions .XRA and .CWN will be
175used to attempt to read the file, producing a number of warning
176messages. To specify that only the GSAS powder data importer be
177used, include a third argument, <tt>fmthint="GSAS powder"</tt> or
178something similar (<tt>fmthint="GSAS"</tt> is also fine, but will
179cause both the "GSAS powder data" and the "GSAS-II .gpx" importer to be considered.)
180<UL><LI>
181  Note that <I>project</I><tt>.add_powder_histogram()</tt> returns a
182  powder histogram objects, which here are saved for later reference. It is
183  also possible to obtain these using <tt>gpx.histograms()</tt>, which
184  returns a list of defined histograms.
185</UL>
186<P>
187Then we add a phase to the project using
188  <I>project</I><tt>.add_phase()</tt>. This specifies a CIF containing the
189  structural information, a name for the phase and specifies that the
190  two histograms are "added" (linked) to the phase. The
191<tt>fmthint="CIF"</tt> parameter can also optionally be specified to limit the
192importers that will be tried.
193
194  <UL><LI>
195  Note that <I>project</I><tt>.add_phase()</tt>
196  returns a phase object</li>
197  <LI>Also, the previously saved histogram objects are used in the
198  <I>project</I><tt>.add_phase()</tt> call. Note that it is
199  also possible to reference the two histgrams by their numbers in the
200  project (here <tt>histograms=[0,1]</tt>) or by the histogram names
201  (here <BR><TT>histograms=['PWDR PBSO4.XRA Bank 1', 'PWDR PBSO4.CWN
202  Bank 1']</TT>).
203  These three ways to use the <tt>histograms</tt> parameter produce
204  the same result.
205</blockquote>
206   
207<A name="ChangeCycles">
208<h4>4: Change the number of refinement cycles</H4>
209</A><blockquote>
210While this is not noted in the original tutorial, this exercise will
211  not complete properly if more variables are added to the refinement
212  without converging the refinement (or at least coming close to
213  convergence) at each refinement step. This is best accomplished by
214  increasing the number of least-squares cycles. No supplied method (at
215  present) allows this parameter to be set straightforwardly, but this
216  code will do this:
217 
218 <blockquote><textarea rows="3" cols="70" readonly>
219# not in tutorial: increase # of cycles to improve convergence
220gpx.data['Controls']['data']['max cyc'] = 8 # not in API </textarea></blockquote>
221
222<P>
223  To find this parameter in the GSAS-II data structure, I followed
224  these steps: In the GUI, Controls is the tree item
225  corresponding to the section where Least Squares cycles are
226  set. Executing the command <tt>gpx.data.keys()</tt> shows the names of the
227  entries in the dictionary corresponding to the GSAS-II tree and
228  indeed one entry is 'Controls'.
229  Command <tt>gpx.data['Controls'].keys()</tt> shows that all
230  values are located in an entry labeled 'data' and
231  <tt>gpx.data['Controls']['data'].keys()</tt> shows the entries in
232this section. Examination of
233<tt>gpx.data['Controls']['data']['max cyc']</tt> shows the value 3,
234which is default number of cycles. Thus,
235the number of cycles is changed with the Python command above.
236
237</blockquote>
238
239<h4>5: Set initial variables and refine</H4>
240<blockquote>
241
242In Step 4 of the original tutorial, the refinement is performed with
243the variables that are on by default. These variables are the two
244histogram scale factors and three background parameters for each
245histogram (8 in total). With GSASIIscriptable, the background
246refinement flags are not turned on by default, so a dictionary must be
247created to set these histogram variables. This code:
248<blockquote><textarea rows="5" cols="75" readonly>
249# tutorial step 4: turn on background refinement (Hist)
250refdict0 = {"set": {"Background": { "no. coeffs": 3, "refine": True }}}
251gpx.save('step4.gpx')
252gpx.do_refinements([refdict0])
253HistStats(gpx)</textarea></blockquote>
254<UL>
255<LI>defines a dictionary (refdict0) that will be used to set the number of background coefficients (not really
256needed, since the default is 3) and sets the background refinement
257flag "on".</li>
258<LI>The project is saved under a new name, so that we can later look at
259the results from each step separately.</li>
260<LI>The parameters in refdict0 are set using
261<I>project</i><tt>.do_refinements()</tt> for both histograms in the
262project and the a refinement is run.
263<LI>The weighted profile R-factor values for all histograms in the project
264are printed using the previously defined <tt>HistStats()</tt> function.
265</UL>
266
267<A name="HistStats">
268Function <tt>HistStats()</tt> works by using <tt>gpx.histograms()</tt>
269to iterate over all histograms in the project, setting <tt>hist</tt> to each histogram
270object. Class member <tt>hist.name</tt> provides the histogram name
271and method <tt>hist.get_wR()</tt> looks up the profile R-factor and
272prints them. The function also writes the final refinement results
273into the current project file. </A>
274
275The output from this will be:<blockquote><pre>
276 Hessian Levenburg-Marquardt SVD refinement on 8 variables:
277initial chi^2 9.6912e+06
278 Cycle: 0, Time: 1.88s, Chi**2: 6.7609e+05, Lambda: 0.001,  Delta: 0.93
279initial chi^2 6.7609e+05
280 Cycle: 1, Time: 1.84s, Chi**2: 6.7602e+05, Lambda: 0.001,  Delta: 0.000104
281converged
282Found 0 SVD zeros
283Read from file: /Users/toby/software/G2/Tutorials/PythonScript/step4.bak0.gpx
284Save to file  : /Users/toby/software/G2/Tutorials/PythonScript/step4.gpx
285GPX file save successful
286 Refinement results are in file: /Users/toby/software/G2/Tutorials/PythonScript/step4.lst
287 ***** Refinement successful *****
288*** profile Rwp, step4.gpx
289        PWDR PBSO4.XRA Bank 1: 40.88
290        PWDR PBSO4.CWN Bank 1: 18.65
291</pre></blockquote>Note that the Rwp values agree with what is
292 expected from the original tutorial.
293<P>
294<B>Note:</B> that there are several equivalent ways to set the histogram
295parameters using
296<tt>G2PwdrData.set_refinements()</tt>,
297<tt>G2Project.set_refinement()</tt> or
298<tt>my_project.do_refinements()</tt>, as
299<A href="http://gsas-ii.readthedocs.io/en/latest/GSASIIscripts.html#refinement-parameters">
300described here</A>. Thus, the <tt>gpx.do_refinements([refdict0])</tt>
301statement above could be replaced with:
302<blockquote><pre>
303gpx.set_refinement({"Background": { "no. coeffs": 3, "refine": True }})
304gpx.do_refinements([{}])
305</pre></blockquote>
306or
307<blockquote><pre> 
308for hist in gpx.histograms():
309    hist.set_refinements({"Background": { "no. coeffs": 3, "refine": True }})
310gpx.do_refinements([{}])
311</pre></blockquote>
312
313</blockquote>
314<h4>6: Add unit cell parameters to refinement</h4>
315<blockquote>
316
317In Step 5 of the original tutorial, the refinement is performed again,
318but the unit cell is refined in the phase.
319
320<blockquote><textarea rows="6" cols="75" readonly>
321# tutorial step 5: add unit cell refinement (Phase)
322gpx.save('step5.gpx')
323refdict1 = {"set": {"Cell": True}} # set the cell flag (for all phases)
324gpx.set_refinement(refdict1)
325gpx.do_refinements([{}])
326HistStats(gpx)</textarea></blockquote>
327
328In this case the <tt>gpx.set_refinement(refdict1)</tt> statement can
329be replaced with:<blockquote><pre>
330phase0.set_refinements({"Cell": True})
331</pre></blockquote>
332
333Note that it is also possible to combine the refinement in the current
334and previous section using
335<blockquote><pre>
336gpx.do_refinements([refdict0,refdict1])
337</pre></blockquote>
338but then the results are not saved in separate project files and it is
339not possible to see the Rwp values after each refinement.
340
341</blockquote>
342<h4>7: Add hydrostatic strain for just one histogram</h4>
343<blockquote>
344
345In Step 6 of the original tutorial, the refinement is performed again
346after adding Histogram/Phase (HAP) parameters so that the lattice
347constants for the first histogram (only) can vary. This is done with
348this code:
349
350<blockquote><textarea rows="6" cols="75" readonly>
351# tutorial step 6: add Dij terms (HAP) for histogram 1 only
352gpx.save('step6.gpx')
353refdict2 = {"set": {"HStrain": True}} # set HAP parameters
354gpx.set_refinement(refdict2,phase=phase0,histogram=[hist1])
355gpx.do_refinements([{}]) # refine after setting
356HistStats(gpx)</textarea></blockquote>
357
358Here we cannot use <tt>gpx.do_refinements([refdict2])</tt> with
359<tt>refdict2</tt> defined as above because
360that would turn on refinement of the Hstrain terms for all histograms
361and all phases. There are several ways to restrict the parameter
362changes to specified histogram(s) and phase(s). One is to call a
363method in the phase object(s) directly, such as
364replacing the <tt>gpx.set_refinement(...)</tt>
365statement with this:
366<blockquote><pre>
367phase0.set_HAP_refinements({"HStrain": True},histograms=[hist1])
368</pre></blockquote>
369
370It is also possible to add "histograms" and/or "phases" values into
371the <tt>refdict2</tt> dict, as will be <a href="#SingeStep">described below.</a>
372</blockquote>
373
374<h4>8: Add X-ray Sample broadening terms</h4>
375<blockquote>
376The next step in the original tutorial is to treat sample broadening
377by turning on refinement of the "Mustrain" (microstrain) and
378"Size" (Scherrer broadening) terms using an isotropic (single-term)
379model. As described in the
380<A
381href="http://gsas-ii.readthedocs.io/en/latest/GSASIIscripts.html#histogram-and-phase-parameters">
382documentation for Histogram/Phase parameters</A>, type must always be
383specfied even as in this case where it is not being changed from the
384existing default. Again, since these parameters are being set only for
385one histogram, either <tt>phase0.set_HAP_refinements()</tt> or
386<tt>gpx.set_refinement()</tt> must be called or to use
387<tt>gpx.do_refinements([refdict3])</tt> the "histograms" element must be
388included inside <tt>refdict3</tt>.
389
390<blockquote><textarea rows="8" cols="75" readonly>
391# tutorial step 7: add size & strain broadening (HAP) for histogram 1 only
392gpx.save('step7.gpx')
393refdict3 = {"set": {"Mustrain": {"type":"isotropic","refine":True},
394                    "Size":{"type":"isotropic","refine":True},
395                    }}
396gpx.set_refinement(refdict3,phase=phase0,histogram=[hist1])
397gpx.do_refinements([{}]) # refine after setting
398HistStats(gpx)</textarea></blockquote>
399
400</blockquote>
401<h4>9: Add Structural and Sample Displacement Parameters</h4>
402<blockquote>
403The step 8 in the original tutorial is to treat sample displacement
404for each histogram/phase. These parameters are different because of
405the differing diffractometer geometries. We also add refinement of
406sample parameters using <tt>phase0.set_refinements()</tt> to set the
407"X" (coordinate) and "U" (displacement) refinement flags for all
408atoms. This is done with this code:
409
410<blockquote><textarea rows="9" cols="75" readonly>
411# tutorial step 8: add sample parameters & set radius (Hist); refine atom parameters (phase)
412gpx.save('step8.gpx')
413hist1.set_refinements({'Sample Parameters': ['Shift']})
414hist2.set_refinements({'Sample Parameters': ['DisplaceX', 'DisplaceY']})
415hist2.data['Sample Parameters']['Gonio. radius'] = 650. # not in API
416phase0.set_refinements({"Atoms":{"all":"XU"}})
417gpx.do_refinements([{}]) # refine after setting
418HistStats(gpx)</textarea></blockquote>
419
420Note that the original tutorial
421calls for the diffractometer radius to be changed to the correct value
422so that the displacement value is in the correct units. This parameter
423cannot be set from any GSASIIscriptable routines, but following a
424similar process,
425<A href="#ChangeCycles">as before</A>, the location for this
426setting can be located in the histogram's 'Sample Parameters' section
427(as <tt>hist2.data['Sample Parameters']['Gonio. radius']</tt>).
428
429Also note that the settings provided in the
430<tt>phase0.set_refinements()</tt> and statements
431<tt>gpx.do_refinements()</tt> could have
432been combined into this single statement:
433<blockquote><pre>
434gpx.do_refinements({"set":{"Atoms":{"all":"XU"}})
435</pre></blockquote>
436
437</blockquote>
438<h4>10: Change Data Limits; Vary Gaussian Profile Terms</h4>
439<blockquote>
440The final step in the original tutorial is to trim the range of data
441used in the refinement to exclude data where no reflections occur and
442where the peaks are cut off at high angle. Also, additional parameters
443are refined here because the supplied instrumental profile parameters
444are not very accurate descriptions for the datasets that are used
445here. It is not possible to refine the Lorentzian x-ray instrumental
446broadening terms, since this broadening is treated as sample
447broadening. The Lorentzian neutron broadening is negligible.
448
449<blockquote><textarea rows="6" cols="75" readonly>
450# tutorial step 9: change data limits & inst. parm refinements (Hist)
451gpx.save('step9.gpx')
452hist1.set_refinements({'Limits': [16.,158.4]})
453hist2.set_refinements({'Limits': [19.,153.]})
454gpx.do_refinements([{"set": {'Instrument Parameters': ['U', 'V', 'W']}}])
455HistStats(gpx)</textarea></blockquote>
456
457Note that the final <tt>gpx.do_refinements()</tt> statement can be
458replaced with calls to
459<tt>hist</tt><I>X</I><tt>.set_refinements()</tt> or
460<tt>gpx.set_refinement()</tt>, such as
461<blockquote><pre>
462hist1.set_refinements({'Instrument Parameters': ['U', 'V', 'W']})
463hist2.set_refinements({'Instrument Parameters': ['U', 'V', 'W']})
464gpx.do_refinements([{}])
465</pre></blockquote>
466
467</blockquote>
468<a name="SingleStep">
469<h2>B. Single-Step Refinement</h2></a>
470<blockquote>
471
472As is noted in the
473<A href="http://gsas-ii.readthedocs.io/en/latest/GSASIIscripts.html" target="_blank">
474GSASIIscriptable module documentation</A>,
475the <I>project</i><tt>.do_refinements()</tt> method can be used to
476perform multiple refinement steps.
477Note that this version of the exercise can be
478<A href="https://subversion.xray.aps.anl.gov/trac/pyGSAS/browser/Tutorials/PythonScript/data/SingleStep.py?format=txt">
479downloaded here</A>.
480To duplicate the above steps into a
481single call, a more complex set of dicts must be created, as shown
482below:
483
484<blockquote><textarea rows="38" cols="75" readonly>
485# tutorial step 4: turn on background refinement (Hist)
486refdict0 = {"set": {"Background": { "no. coeffs": 3, "refine": True }},
487            "output":'step4.gpx',
488            "call":HistStats,"callargs":[gpx]} # callargs actually unneeded
489                                               # as [gpx] is the default
490# tutorial step 5: add unit cell refinement (Phase)
491refdict1 = {"set": {"Cell": True},    # set the cell flag (for all phases)
492            "output":'step5.gpx', "call":HistStats}
493# tutorial step 6: add Dij terms (HAP) for phase 1 only
494refdict2 = {"set": {"HStrain": True},  # set HAP parameters
495            "histograms":[hist1],      # histogram 1 only
496            "phases":[phase0],         # unneeded (default is all
497                                       # phases) included as a example
498            "output":'step6.gpx', "call":HistStats}
499# tutorial step 7: add size & strain broadening (HAP) for histogram 1 only
500refdict3 = {"set": {"Mustrain": {"type":"isotropic","refine":True},
501                    "Size":{"type":"isotropic","refine":True},},
502            "histograms":[hist1],      # histogram 1 only
503            "output":'step7.gpx', "call":HistStats}
504# tutorial step 8: add sample parameters & set radius (Hist); refine
505#                  atom parameters (phase)
506refdict4a = {"set": {'Sample Parameters': ['Shift']},
507            "histograms":[hist1],      # histogram 1 only
508             "skip": True}
509refdict4b = {"set": {"Atoms":{"all":"XU"}, # not affected by histograms
510             'Sample Parameters': ['DisplaceX', 'DisplaceY']},
511            "histograms":[hist2],      # histogram 2 only
512            "output":'step8.gpx', "call":HistStats}
513# tutorial step 9: change data limits & inst. parm refinements (Hist)
514refdict5a = {"set": {'Limits': [16.,158.4]},
515            "histograms":[hist1],      # histogram 1 only
516             "skip": True,}
517refdict5b = {"set": {'Limits': [19.,153.]},
518            "histograms":[hist2],      # histogram 2 only
519             "skip": True}
520refdict5c = {"set": {'Instrument Parameters': ['U', 'V', 'W']},
521            "output":'step9.gpx', "call":HistStats}
522</textarea></blockquote>
523
524Note that above the "call" and "callargs" dict entries are
525defined to run <tt>HistStats</tt> and "output" is used to designate
526the .gpx file that will be needed. When parameters should be changed
527in specific histograms, the entries <tt>"histograms":[hist1]</tt> and
528<tt>"histograms":[hist2]</tt> are used
529(equivalent would be <tt>"histograms":[0]</tt> and
530<tt>"histograms":[1]</tt>).
531Since there is only one phase present, use of <tt>"phase":[0]</tt> is
532superfluous, but in a more complex refinement, this could be needed.
533<p>
534A list of dicts is then prepared here:
535<blockquote><textarea rows="3" cols="75" readonly>
536dictList = [refdict0,refdict1,refdict2,refdict3,
537            refdict4a,refdict4b,
538            refdict5a,refdict5b,refdict5c]</textarea></blockquote>
539
540Steps 4 through 10, above then can be performed with these few commands:
541<blockquote><textarea rows="4" cols="75" readonly>
542# Change number of cycles and radius
543gpx.data['Controls']['data']['max cyc'] = 8 # not in API
544hist2.data['Sample Parameters']['Gonio. radius'] = 650. # not in API
545gpx.do_refinements(dictList)</textarea></blockquote>
546
547</blockquote>
548<a name="Simulate">
549<h2>C. Powder Pattern Simulation</h2></a>
550<blockquote>
551Use of the
552<A href="http://gsas-ii.readthedocs.io/en/latest/GSASIIscripts.html" target="_blank">
553GSASIIscriptable module</A> makes it very simple to simulate a powder
554diffraction pattern using GSAS-II. This is demonstrated in this
555section.
556
557Note that the Python commands can be
558<A
559href="https://subversion.xray.aps.anl.gov/trac/pyGSAS/browser/Tutorials/PythonScript/data/sim.py?format=txt">downloaded here</A>.
560<P>
561As before, the location of the GSASIIscripts Python module must be
562defined and the module must be loaded:
563
564<blockquote><textarea rows="3" cols="60" readonly>
565import os,sys
566sys.path.insert(0,os.path.expanduser("~/g2conda/GSASII/"))
567import GSASIIscriptable as G2sc</textarea></blockquote>
568
569To simplify this example, as before we will define the location where files will
570be read from and written (as <tt>datadir</tt> and
571<tt>workdir</tt>). Note that files "inst_d1a.prm" and
572"PbSO4-Wyckoff.cif" from <A
573href="https://subversion.xray.aps.anl.gov/pyGSAS/Tutorials/PythonScript/data/">here</A>
574are needed.
575
576<blockquote><textarea rows="2" cols="70" readonly>
577workdir = "/Users/toby/Scratch/PythonScript"
578datadir = "/Users/toby/software/G2/Tutorials/PythonScript/data"</textarea></blockquote>
579
580We then need to create a project and for this example we choose to
581define the phase first. (It would work equally well to create the
582histogram first and then define the phase.)
583
584<blockquote><textarea rows="4" cols="70" readonly>
585gpx = G2sc.G2Project(newgpx='PbSO4sim.gpx') # create a project
586# step 1, setup: add a phase to the project
587phase0 = gpx.add_phase(os.path.join(datadir,"PbSO4-Wyckoff.cif"),
588                      phasename="PbSO4",fmthint='CIF') </textarea></blockquote>
589
590We then add a "dummy" histogram to the project. Note that an
591instrument parameter file is specified, but not a data file. The range
592of data to be used and the step size must be specified. The phases
593parameter is specified as ``gpx.phases()`` which creates a list of all
594the previously read phases, which in this case is equivalent to
595``[phase0]``.
596
597<blockquote><textarea rows="6" cols="80" readonly>
598# step 2, setup: add a simulated histogram and link it to the previous phase(s)
599hist1 = gpx.add_simulated_powder_histogram("PbSO4 simulation",
600                          os.path.join(datadir,"inst_d1a.prm"),
601                          5.,120.,0.01,
602                          phases=gpx.phases())</textarea></blockquote>
603
604Note that the computed pattern cannot be seen above "simulated noise"
605unless the intensities are large enough. We can change the pattern
606scale factor using the Scale factor (parameter
607<tt>hist1.data['Sample Parameters']['Scale'][0]</tt>), as shown below.
608
609<blockquote><textarea rows="3" cols="70" readonly>
610# Step 3: Set the scale factor to adjust the y scale
611hist1.SampleParameters['Scale'][0] = 1000000.
612</textarea></blockquote>
613
614Next, to perform the simulation computation, a refinement is
615needed:
616
617<blockquote><textarea rows="4" cols="80" readonly>
618# step 4, compute: turn off parameter optimization and calculate pattern
619gpx.data['Controls']['data']['max cyc'] = 0 # refinement not needed
620gpx.do_refinements([{}])
621gpx.save()</textarea></blockquote>
622
623However, there is no need to actually optimize any variables,
624so the number of refinement cycles is set to zero. Refinement is
625initiated then with <I>proj</i>.<tt>do_refinements</tt>. To keep the
626results in a .gpx file, the
627project is saved.
628
629<P>
630Finally, we want to do something with the results. The histogram could
631be written to a file with the <A href=
632"http://gsas-ii.readthedocs.io/en/latest/GSASIIscriptable.html#GSASIIscriptable.G2PwdrData.Export">
633<i>histogram.</i><tt>Export()</tt></A> command. Note that the first time
634a refinement computation is done with a dummy histogram the "observed"
635pattern is set from the calculated intensities. Here an alternate option is
636used, where the values are retrieved and plotted.
637
638<blockquote><textarea rows="6" cols="70" readonly>
639# step 5, retrieve results & plot
640x = gpx.histogram(0).getdata('x')
641y = gpx.histogram(0).getdata('ycalc')
642import matplotlib.pyplot as plt
643plt.plot(x,y)
644plt.savefig('PbSO4.png') # to show on screen use: plt.show()</textarea></blockquote>
645
646Note that in the above, <tt>gpx.histogram(0)</tt> is used to show how
647to reference <tt>hist1</tt> when the histogram object is not
648known. The generated two-theta values and computed intensity values
649are retrieved and the remaining lines generate a very simple plot
650which is saved to a file, as shown below.
651</blockquote>
652
653<img src="PbSO4.png">
654
655<hr>
656<address></address>
657<!-- hhmts start -->Last modified: Mon Feb 19 14:53:43 CST 2018 <!-- hhmts end -->
658</body> </html>
Note: See TracBrowser for help on using the repository browser.