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> |
---|
8 | To demonstrate the use of the |
---|
9 | <A href="http://gsas-ii.readthedocs.io/en/latest/GSASIIscripts.html" target="_blank"> |
---|
10 | GSASIIscriptable module</A>. This uses a Python script |
---|
11 | to perform a refinement or computation, but without use of the GSAS-II |
---|
12 | graphical user interface. Note that the |
---|
13 | <A href="http://gsas-ii.readthedocs.io/en/latest/GSASIIscripts.html" target="_blank"> |
---|
14 | GSASIIscriptable module</A> does not offer access to all of the |
---|
15 | features within GSAS-II, but over time it is expected to grow in |
---|
16 | capabilities. Much of the initial development for this module was done |
---|
17 | by Jackson O'Donnell, as a summer undergraduate visitor under |
---|
18 | supervisor Dr. Maria Chan. Other programming contributions are welcome. |
---|
19 | <P> |
---|
20 | This 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> |
---|
26 | The first section duplicates the refinement in the |
---|
27 | <A |
---|
28 | href="https://subversion.xray.aps.anl.gov/pyGSAS/Tutorials/CWCombined/Combined%20refinement.htm" |
---|
29 | target="_blank">GSAS-II CW Combined Refinement</A> |
---|
30 | tutorial as a multi-step process. The second section repeats the previous |
---|
31 | refinement, but demonstrates how a complex process can be entered into |
---|
32 | a single Python dict and the refinement executed with a single |
---|
33 | function call. The third section shows how a pattern simulation, |
---|
34 | rather than refinement can be executed from a script. |
---|
35 | |
---|
36 | |
---|
37 | <h2>Prerequisites</h2> |
---|
38 | This exercise assumes that the reader has reasonable familiarity with |
---|
39 | the Python language. |
---|
40 | Also, the documentation on the |
---|
41 | <A href="http://gsas-ii.readthedocs.io/en/latest/GSASIIscripts.html" |
---|
42 | target="_blank"> |
---|
43 | Scripting Tools in the GSAS-II GSASIIscriptable module</A> should be |
---|
44 | read from here: <A href="http://gsas-ii.readthedocs.io/en/latest/GSASIIscripts.html" |
---|
45 | target="_blank">http://gsas-ii.readthedocs.io/en/latest/GSASIIscripts.html</A>. It |
---|
46 | is also wise to read the <A href="https://subversion.xray.aps.anl.gov/pyGSAS/Tutorials/CWCombined/Combined%20refinement.htm" |
---|
47 | target="_blank">GSAS-II CW Combined Refinement</A> tutorial that this |
---|
48 | exercise is modeled upon, which explains why each refinement step is |
---|
49 | being used. |
---|
50 | The exercise will require downloading of the same files needed for the |
---|
51 | Combined Refinement tutorial: "PBSO4.XRA", "INST_XRY.PRM", "PBSO4.CWN", |
---|
52 | "inst_d1a.prm" and "PbSO4-Wyckoff.cif", |
---|
53 | which can be downloaded from |
---|
54 | <A href="https://subversion.xray.aps.anl.gov/pyGSAS/Tutorials/PythonScript/data/"> |
---|
55 | https://subversion.xray.aps.anl.gov/pyGSAS/Tutorials/PythonScript/data/</A>) |
---|
56 | or will be downloaded with this tutorial if that is requested. |
---|
57 | <P> |
---|
58 | The exercise can be performed by placing all of the Python commands |
---|
59 | into a script, |
---|
60 | but a more pedagogical approach will be to enter the |
---|
61 | commands into a Python interpreter. Use of IPython or Jupyter to run |
---|
62 | Python will make this a more pleasant experience. |
---|
63 | |
---|
64 | <a name="multistep"> |
---|
65 | <h2>A. Create a Multi-step Script</h2></a> |
---|
66 | |
---|
67 | Note 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"> |
---|
69 | download here</A>, but some editing will be needed to match where files |
---|
70 | are found on your computer. |
---|
71 | |
---|
72 | <h4>0: Load the GSASIIscriptable module</H4> |
---|
73 | <blockquote> |
---|
74 | The first step in writing any Python module is to load the modules |
---|
75 | that will be needed. Here we will need the <tt>os</tt> and |
---|
76 | <tt>sys</tt> modules from the Python standard library and the |
---|
77 | GSASIIscriptable from the GSAS-II code. The location where the |
---|
78 | GSAS-II source code is installed must usually be hard-coded into the |
---|
79 | script as is done in the example below. Note that a common location |
---|
80 | for this will be |
---|
81 | <TT>os.path.expanduser("~/g2conda/GSASII/")</TT> or |
---|
82 | <tt>'/Users/toby/software/GSASII'</tt>, etc. |
---|
83 | Thus the script will begin with something like this: |
---|
84 | |
---|
85 | <blockquote><textarea rows="3" cols="60" readonly> |
---|
86 | import os,sys |
---|
87 | sys.path.insert(0,os.path.expanduser("~/g2conda/GSASII/")) |
---|
88 | import GSASIIscriptable as G2sc</textarea></blockquote> |
---|
89 | |
---|
90 | If a <tt>ImportError: No module named GSASIIscriptable</tt> error |
---|
91 | occurs, then the wrong directory location has been supplied for the |
---|
92 | GSAS-II files. |
---|
93 | </blockquote> |
---|
94 | |
---|
95 | <h4>1: Define some other prerequisite code</H4> |
---|
96 | <blockquote> |
---|
97 | To simplify this example, we will define the location where files will |
---|
98 | be written as <tt>workdir</tt> (this directory must exist, but it may |
---|
99 | be 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 |
---|
102 | be found as <tt>datadir</tt>. (As discussed previously, |
---|
103 | these files can be |
---|
104 | <A href="https://subversion.xray.aps.anl.gov/pyGSAS/Tutorials/PythonScript/data/"> |
---|
105 | downloaded from here. </A>) |
---|
106 | <P> |
---|
107 | We also define here a short function to display the weighted profile R |
---|
108 | factor for every histogram in a project, <tt>HistStats</tt>. This will |
---|
109 | be <A href="#HistStats">discussed later</A>. |
---|
110 | |
---|
111 | <blockquote><textarea rows="10" cols="70" readonly> |
---|
112 | workdir = "/Users/toby/Scratch/PythonScript" |
---|
113 | datadir = "/Users/toby/software/G2/Tutorials/PythonScript/data" |
---|
114 | |
---|
115 | def 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> |
---|
127 | The first step in creating a GSASIIscriptable script to to create or |
---|
128 | access a GSAS-II project, which is done with a call to |
---|
129 | <TT>GSASIIscriptable.G2Project()</tt>. This can be done in one of two |
---|
130 | ways: 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 |
---|
133 | reads an existing project (.gpx) file. Both return a <tt>G2Project</tt> wrapper |
---|
134 | object that is used to access a number of methods and variables. |
---|
135 | Note that <TT>GSASIIscriptable.G2Project()</tt> can read .gpx files |
---|
136 | written by both previous <TT>GSASIIscriptable.G2Project()</tt> runs or |
---|
137 | the GSAS-II GUI. |
---|
138 | <P> |
---|
139 | In 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 |
---|
144 | gpx = G2sc.G2Project(newgpx='PbSO4.gpx')</textarea></blockquote> |
---|
145 | |
---|
146 | Note that the file will not actually be created until an operation |
---|
147 | that saves the project is called. |
---|
148 | </blockquote> |
---|
149 | |
---|
150 | <h4>3: Add Histograms and Phase to the GSAS-II project</H4> |
---|
151 | <blockquote> |
---|
152 | To add two powder diffraction datasets (histograms) and a phase to the |
---|
153 | project using this code: |
---|
154 | |
---|
155 | <blockquote><textarea rows="10" cols="70" readonly> |
---|
156 | # setup step 1: add two histograms to the project |
---|
157 | hist1 = gpx.add_powder_histogram(os.path.join(datadir,"PBSO4.XRA"), |
---|
158 | os.path.join(datadir,"INST_XRY.PRM")) |
---|
159 | hist2 = 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 |
---|
162 | phase0 = gpx.add_phase(os.path.join(datadir,"PbSO4-Wyckoff.cif"), |
---|
163 | phasename="PbSO4", |
---|
164 | histograms=[hist1,hist2])</textarea></blockquote> |
---|
165 | |
---|
166 | |
---|
167 | We use the <I>project</I><tt>.add_powder_histogram()</tt> method in |
---|
168 | the <TT>GSASIIscriptable</tt> class to read in the powder data. |
---|
169 | The 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. |
---|
174 | All importers that allow for extensions .XRA and .CWN will be |
---|
175 | used to attempt to read the file, producing a number of warning |
---|
176 | messages. To specify that only the GSAS powder data importer be |
---|
177 | used, include a third argument, <tt>fmthint="GSAS powder"</tt> or |
---|
178 | something similar (<tt>fmthint="GSAS"</tt> is also fine, but will |
---|
179 | cause 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> |
---|
187 | Then 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 |
---|
192 | importers 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> |
---|
210 | While 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 |
---|
220 | gpx.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 |
---|
232 | this section. Examination of |
---|
233 | <tt>gpx.data['Controls']['data']['max cyc']</tt> shows the value 3, |
---|
234 | which is default number of cycles. Thus, |
---|
235 | the 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 | |
---|
242 | In Step 4 of the original tutorial, the refinement is performed with |
---|
243 | the variables that are on by default. These variables are the two |
---|
244 | histogram scale factors and three background parameters for each |
---|
245 | histogram (8 in total). With GSASIIscriptable, the background |
---|
246 | refinement flags are not turned on by default, so a dictionary must be |
---|
247 | created to set these histogram variables. This code: |
---|
248 | <blockquote><textarea rows="5" cols="75" readonly> |
---|
249 | # tutorial step 4: turn on background refinement (Hist) |
---|
250 | refdict0 = {"set": {"Background": { "no. coeffs": 3, "refine": True }}} |
---|
251 | gpx.save('step4.gpx') |
---|
252 | gpx.do_refinements([refdict0]) |
---|
253 | HistStats(gpx)</textarea></blockquote> |
---|
254 | <UL> |
---|
255 | <LI>defines a dictionary (refdict0) that will be used to set the number of background coefficients (not really |
---|
256 | needed, since the default is 3) and sets the background refinement |
---|
257 | flag "on".</li> |
---|
258 | <LI>The project is saved under a new name, so that we can later look at |
---|
259 | the 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 |
---|
262 | project and the a refinement is run. |
---|
263 | <LI>The weighted profile R-factor values for all histograms in the project |
---|
264 | are printed using the previously defined <tt>HistStats()</tt> function. |
---|
265 | </UL> |
---|
266 | |
---|
267 | <A name="HistStats"> |
---|
268 | Function <tt>HistStats()</tt> works by using <tt>gpx.histograms()</tt> |
---|
269 | to iterate over all histograms in the project, setting <tt>hist</tt> to each histogram |
---|
270 | object. Class member <tt>hist.name</tt> provides the histogram name |
---|
271 | and method <tt>hist.get_wR()</tt> looks up the profile R-factor and |
---|
272 | prints them. The function also writes the final refinement results |
---|
273 | into the current project file. </A> |
---|
274 | |
---|
275 | The output from this will be:<blockquote><pre> |
---|
276 | Hessian Levenburg-Marquardt SVD refinement on 8 variables: |
---|
277 | initial chi^2 9.6912e+06 |
---|
278 | Cycle: 0, Time: 1.88s, Chi**2: 6.7609e+05, Lambda: 0.001, Delta: 0.93 |
---|
279 | initial chi^2 6.7609e+05 |
---|
280 | Cycle: 1, Time: 1.84s, Chi**2: 6.7602e+05, Lambda: 0.001, Delta: 0.000104 |
---|
281 | converged |
---|
282 | Found 0 SVD zeros |
---|
283 | Read from file: /Users/toby/software/G2/Tutorials/PythonScript/step4.bak0.gpx |
---|
284 | Save to file : /Users/toby/software/G2/Tutorials/PythonScript/step4.gpx |
---|
285 | GPX 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 |
---|
295 | parameters 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"> |
---|
300 | described here</A>. Thus, the <tt>gpx.do_refinements([refdict0])</tt> |
---|
301 | statement above could be replaced with: |
---|
302 | <blockquote><pre> |
---|
303 | gpx.set_refinement({"Background": { "no. coeffs": 3, "refine": True }}) |
---|
304 | gpx.do_refinements([{}]) |
---|
305 | </pre></blockquote> |
---|
306 | or |
---|
307 | <blockquote><pre> |
---|
308 | for hist in gpx.histograms(): |
---|
309 | hist.set_refinements({"Background": { "no. coeffs": 3, "refine": True }}) |
---|
310 | gpx.do_refinements([{}]) |
---|
311 | </pre></blockquote> |
---|
312 | |
---|
313 | </blockquote> |
---|
314 | <h4>6: Add unit cell parameters to refinement</h4> |
---|
315 | <blockquote> |
---|
316 | |
---|
317 | In Step 5 of the original tutorial, the refinement is performed again, |
---|
318 | but 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) |
---|
322 | gpx.save('step5.gpx') |
---|
323 | refdict1 = {"set": {"Cell": True}} # set the cell flag (for all phases) |
---|
324 | gpx.set_refinement(refdict1) |
---|
325 | gpx.do_refinements([{}]) |
---|
326 | HistStats(gpx)</textarea></blockquote> |
---|
327 | |
---|
328 | In this case the <tt>gpx.set_refinement(refdict1)</tt> statement can |
---|
329 | be replaced with:<blockquote><pre> |
---|
330 | phase0.set_refinements({"Cell": True}) |
---|
331 | </pre></blockquote> |
---|
332 | |
---|
333 | Note that it is also possible to combine the refinement in the current |
---|
334 | and previous section using |
---|
335 | <blockquote><pre> |
---|
336 | gpx.do_refinements([refdict0,refdict1]) |
---|
337 | </pre></blockquote> |
---|
338 | but then the results are not saved in separate project files and it is |
---|
339 | not 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 | |
---|
345 | In Step 6 of the original tutorial, the refinement is performed again |
---|
346 | after adding Histogram/Phase (HAP) parameters so that the lattice |
---|
347 | constants for the first histogram (only) can vary. This is done with |
---|
348 | this code: |
---|
349 | |
---|
350 | <blockquote><textarea rows="6" cols="75" readonly> |
---|
351 | # tutorial step 6: add Dij terms (HAP) for histogram 1 only |
---|
352 | gpx.save('step6.gpx') |
---|
353 | refdict2 = {"set": {"HStrain": True}} # set HAP parameters |
---|
354 | gpx.set_refinement(refdict2,phase=phase0,histogram=[hist1]) |
---|
355 | gpx.do_refinements([{}]) # refine after setting |
---|
356 | HistStats(gpx)</textarea></blockquote> |
---|
357 | |
---|
358 | Here we cannot use <tt>gpx.do_refinements([refdict2])</tt> with |
---|
359 | <tt>refdict2</tt> defined as above because |
---|
360 | that would turn on refinement of the Hstrain terms for all histograms |
---|
361 | and all phases. There are several ways to restrict the parameter |
---|
362 | changes to specified histogram(s) and phase(s). One is to call a |
---|
363 | method in the phase object(s) directly, such as |
---|
364 | replacing the <tt>gpx.set_refinement(...)</tt> |
---|
365 | statement with this: |
---|
366 | <blockquote><pre> |
---|
367 | phase0.set_HAP_refinements({"HStrain": True},histograms=[hist1]) |
---|
368 | </pre></blockquote> |
---|
369 | |
---|
370 | It is also possible to add "histograms" and/or "phases" values into |
---|
371 | the <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> |
---|
376 | The next step in the original tutorial is to treat sample broadening |
---|
377 | by turning on refinement of the "Mustrain" (microstrain) and |
---|
378 | "Size" (Scherrer broadening) terms using an isotropic (single-term) |
---|
379 | model. As described in the |
---|
380 | <A |
---|
381 | href="http://gsas-ii.readthedocs.io/en/latest/GSASIIscripts.html#histogram-and-phase-parameters"> |
---|
382 | documentation for Histogram/Phase parameters</A>, type must always be |
---|
383 | specfied even as in this case where it is not being changed from the |
---|
384 | existing default. Again, since these parameters are being set only for |
---|
385 | one 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 |
---|
388 | included 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 |
---|
392 | gpx.save('step7.gpx') |
---|
393 | refdict3 = {"set": {"Mustrain": {"type":"isotropic","refine":True}, |
---|
394 | "Size":{"type":"isotropic","refine":True}, |
---|
395 | }} |
---|
396 | gpx.set_refinement(refdict3,phase=phase0,histogram=[hist1]) |
---|
397 | gpx.do_refinements([{}]) # refine after setting |
---|
398 | HistStats(gpx)</textarea></blockquote> |
---|
399 | |
---|
400 | </blockquote> |
---|
401 | <h4>9: Add Structural and Sample Displacement Parameters</h4> |
---|
402 | <blockquote> |
---|
403 | The step 8 in the original tutorial is to treat sample displacement |
---|
404 | for each histogram/phase. These parameters are different because of |
---|
405 | the differing diffractometer geometries. We also add refinement of |
---|
406 | sample parameters using <tt>phase0.set_refinements()</tt> to set the |
---|
407 | "X" (coordinate) and "U" (displacement) refinement flags for all |
---|
408 | atoms. 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) |
---|
412 | gpx.save('step8.gpx') |
---|
413 | hist1.set_refinements({'Sample Parameters': ['Shift']}) |
---|
414 | hist2.set_refinements({'Sample Parameters': ['DisplaceX', 'DisplaceY']}) |
---|
415 | hist2.data['Sample Parameters']['Gonio. radius'] = 650. # not in API |
---|
416 | phase0.set_refinements({"Atoms":{"all":"XU"}}) |
---|
417 | gpx.do_refinements([{}]) # refine after setting |
---|
418 | HistStats(gpx)</textarea></blockquote> |
---|
419 | |
---|
420 | Note that the original tutorial |
---|
421 | calls for the diffractometer radius to be changed to the correct value |
---|
422 | so that the displacement value is in the correct units. This parameter |
---|
423 | cannot be set from any GSASIIscriptable routines, but following a |
---|
424 | similar process, |
---|
425 | <A href="#ChangeCycles">as before</A>, the location for this |
---|
426 | setting can be located in the histogram's 'Sample Parameters' section |
---|
427 | (as <tt>hist2.data['Sample Parameters']['Gonio. radius']</tt>). |
---|
428 | |
---|
429 | Also note that the settings provided in the |
---|
430 | <tt>phase0.set_refinements()</tt> and statements |
---|
431 | <tt>gpx.do_refinements()</tt> could have |
---|
432 | been combined into this single statement: |
---|
433 | <blockquote><pre> |
---|
434 | gpx.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> |
---|
440 | The final step in the original tutorial is to trim the range of data |
---|
441 | used in the refinement to exclude data where no reflections occur and |
---|
442 | where the peaks are cut off at high angle. Also, additional parameters |
---|
443 | are refined here because the supplied instrumental profile parameters |
---|
444 | are not very accurate descriptions for the datasets that are used |
---|
445 | here. It is not possible to refine the Lorentzian x-ray instrumental |
---|
446 | broadening terms, since this broadening is treated as sample |
---|
447 | broadening. 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) |
---|
451 | gpx.save('step9.gpx') |
---|
452 | hist1.set_refinements({'Limits': [16.,158.4]}) |
---|
453 | hist2.set_refinements({'Limits': [19.,153.]}) |
---|
454 | gpx.do_refinements([{"set": {'Instrument Parameters': ['U', 'V', 'W']}}]) |
---|
455 | HistStats(gpx)</textarea></blockquote> |
---|
456 | |
---|
457 | Note that the final <tt>gpx.do_refinements()</tt> statement can be |
---|
458 | replaced 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> |
---|
462 | hist1.set_refinements({'Instrument Parameters': ['U', 'V', 'W']}) |
---|
463 | hist2.set_refinements({'Instrument Parameters': ['U', 'V', 'W']}) |
---|
464 | gpx.do_refinements([{}]) |
---|
465 | </pre></blockquote> |
---|
466 | |
---|
467 | </blockquote> |
---|
468 | <a name="SingleStep"> |
---|
469 | <h2>B. Single-Step Refinement</h2></a> |
---|
470 | <blockquote> |
---|
471 | |
---|
472 | As is noted in the |
---|
473 | <A href="http://gsas-ii.readthedocs.io/en/latest/GSASIIscripts.html" target="_blank"> |
---|
474 | GSASIIscriptable module documentation</A>, |
---|
475 | the <I>project</i><tt>.do_refinements()</tt> method can be used to |
---|
476 | perform multiple refinement steps. |
---|
477 | Note 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"> |
---|
479 | downloaded here</A>. |
---|
480 | To duplicate the above steps into a |
---|
481 | single call, a more complex set of dicts must be created, as shown |
---|
482 | below: |
---|
483 | |
---|
484 | <blockquote><textarea rows="38" cols="75" readonly> |
---|
485 | # tutorial step 4: turn on background refinement (Hist) |
---|
486 | refdict0 = {"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) |
---|
491 | refdict1 = {"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 |
---|
494 | refdict2 = {"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 |
---|
500 | refdict3 = {"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) |
---|
506 | refdict4a = {"set": {'Sample Parameters': ['Shift']}, |
---|
507 | "histograms":[hist1], # histogram 1 only |
---|
508 | "skip": True} |
---|
509 | refdict4b = {"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) |
---|
514 | refdict5a = {"set": {'Limits': [16.,158.4]}, |
---|
515 | "histograms":[hist1], # histogram 1 only |
---|
516 | "skip": True,} |
---|
517 | refdict5b = {"set": {'Limits': [19.,153.]}, |
---|
518 | "histograms":[hist2], # histogram 2 only |
---|
519 | "skip": True} |
---|
520 | refdict5c = {"set": {'Instrument Parameters': ['U', 'V', 'W']}, |
---|
521 | "output":'step9.gpx', "call":HistStats} |
---|
522 | </textarea></blockquote> |
---|
523 | |
---|
524 | Note that above the "call" and "callargs" dict entries are |
---|
525 | defined to run <tt>HistStats</tt> and "output" is used to designate |
---|
526 | the .gpx file that will be needed. When parameters should be changed |
---|
527 | in 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>). |
---|
531 | Since there is only one phase present, use of <tt>"phase":[0]</tt> is |
---|
532 | superfluous, but in a more complex refinement, this could be needed. |
---|
533 | <p> |
---|
534 | A list of dicts is then prepared here: |
---|
535 | <blockquote><textarea rows="3" cols="75" readonly> |
---|
536 | dictList = [refdict0,refdict1,refdict2,refdict3, |
---|
537 | refdict4a,refdict4b, |
---|
538 | refdict5a,refdict5b,refdict5c]</textarea></blockquote> |
---|
539 | |
---|
540 | Steps 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 |
---|
543 | gpx.data['Controls']['data']['max cyc'] = 8 # not in API |
---|
544 | hist2.data['Sample Parameters']['Gonio. radius'] = 650. # not in API |
---|
545 | gpx.do_refinements(dictList)</textarea></blockquote> |
---|
546 | |
---|
547 | </blockquote> |
---|
548 | <a name="Simulate"> |
---|
549 | <h2>C. Powder Pattern Simulation</h2></a> |
---|
550 | <blockquote> |
---|
551 | Use of the |
---|
552 | <A href="http://gsas-ii.readthedocs.io/en/latest/GSASIIscripts.html" target="_blank"> |
---|
553 | GSASIIscriptable module</A> makes it very simple to simulate a powder |
---|
554 | diffraction pattern using GSAS-II. This is demonstrated in this |
---|
555 | section. |
---|
556 | |
---|
557 | Note that the Python commands can be |
---|
558 | <A |
---|
559 | href="https://subversion.xray.aps.anl.gov/trac/pyGSAS/browser/Tutorials/PythonScript/data/sim.py?format=txt">downloaded here</A>. |
---|
560 | <P> |
---|
561 | As before, the location of the GSASIIscripts Python module must be |
---|
562 | defined and the module must be loaded: |
---|
563 | |
---|
564 | <blockquote><textarea rows="3" cols="60" readonly> |
---|
565 | import os,sys |
---|
566 | sys.path.insert(0,os.path.expanduser("~/g2conda/GSASII/")) |
---|
567 | import GSASIIscriptable as G2sc</textarea></blockquote> |
---|
568 | |
---|
569 | To simplify this example, as before we will define the location where files will |
---|
570 | be 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 |
---|
573 | href="https://subversion.xray.aps.anl.gov/pyGSAS/Tutorials/PythonScript/data/">here</A> |
---|
574 | are needed. |
---|
575 | |
---|
576 | <blockquote><textarea rows="2" cols="70" readonly> |
---|
577 | workdir = "/Users/toby/Scratch/PythonScript" |
---|
578 | datadir = "/Users/toby/software/G2/Tutorials/PythonScript/data"</textarea></blockquote> |
---|
579 | |
---|
580 | We then need to create a project and for this example we choose to |
---|
581 | define the phase first. (It would work equally well to create the |
---|
582 | histogram first and then define the phase.) |
---|
583 | |
---|
584 | <blockquote><textarea rows="4" cols="70" readonly> |
---|
585 | gpx = G2sc.G2Project(newgpx='PbSO4sim.gpx') # create a project |
---|
586 | # step 1, setup: add a phase to the project |
---|
587 | phase0 = gpx.add_phase(os.path.join(datadir,"PbSO4-Wyckoff.cif"), |
---|
588 | phasename="PbSO4",fmthint='CIF') </textarea></blockquote> |
---|
589 | |
---|
590 | We then add a "dummy" histogram to the project. Note that an |
---|
591 | instrument parameter file is specified, but not a data file. The range |
---|
592 | of data to be used and the step size must be specified. The phases |
---|
593 | parameter is specified as ``gpx.phases()`` which creates a list of all |
---|
594 | the 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) |
---|
599 | hist1 = 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 | |
---|
604 | Note that the computed pattern cannot be seen above "simulated noise" |
---|
605 | unless the intensities are large enough. We can change the pattern |
---|
606 | scale 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 |
---|
611 | hist1.SampleParameters['Scale'][0] = 1000000. |
---|
612 | </textarea></blockquote> |
---|
613 | |
---|
614 | Next, to perform the simulation computation, a refinement is |
---|
615 | needed: |
---|
616 | |
---|
617 | <blockquote><textarea rows="4" cols="80" readonly> |
---|
618 | # step 4, compute: turn off parameter optimization and calculate pattern |
---|
619 | gpx.data['Controls']['data']['max cyc'] = 0 # refinement not needed |
---|
620 | gpx.do_refinements([{}]) |
---|
621 | gpx.save()</textarea></blockquote> |
---|
622 | |
---|
623 | However, there is no need to actually optimize any variables, |
---|
624 | so the number of refinement cycles is set to zero. Refinement is |
---|
625 | initiated then with <I>proj</i>.<tt>do_refinements</tt>. To keep the |
---|
626 | results in a .gpx file, the |
---|
627 | project is saved. |
---|
628 | |
---|
629 | <P> |
---|
630 | Finally, we want to do something with the results. The histogram could |
---|
631 | be 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 |
---|
634 | a refinement computation is done with a dummy histogram the "observed" |
---|
635 | pattern is set from the calculated intensities. Here an alternate option is |
---|
636 | used, where the values are retrieved and plotted. |
---|
637 | |
---|
638 | <blockquote><textarea rows="6" cols="70" readonly> |
---|
639 | # step 5, retrieve results & plot |
---|
640 | x = gpx.histogram(0).getdata('x') |
---|
641 | y = gpx.histogram(0).getdata('ycalc') |
---|
642 | import matplotlib.pyplot as plt |
---|
643 | plt.plot(x,y) |
---|
644 | plt.savefig('PbSO4.png') # to show on screen use: plt.show()</textarea></blockquote> |
---|
645 | |
---|
646 | Note that in the above, <tt>gpx.histogram(0)</tt> is used to show how |
---|
647 | to reference <tt>hist1</tt> when the histogram object is not |
---|
648 | known. The generated two-theta values and computed intensity values |
---|
649 | are retrieved and the remaining lines generate a very simple plot |
---|
650 | which 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> |
---|