Running a GSAS-II Refinement From a Python Script

In this training example we create a Python script to duplicate the refinement in the GSAS-II CW Combined Refinement tutorial. This uses the GSASIIscriptable module to perform the same refinements steps, but without use of the GSAS-II graphical user interface.

Prerequisites

This exercise assumes that the reader has reasonable familiarity with the Python language. Before beginning this exercise, the documentation on the Scripting Tools in the GSAS-II GSASIIscriptable module should be read from here: http://gsas-ii.readthedocs.io/en/latest/GSASIIscripts.html. It is also wise to read the GSAS-II CW Combined Refinement tutorial that this exercise is modeled upon, which explains why each refinement step is being used.

The exercise can be performed by placing all of the Python commands into a script, but more pedagogical approach will be to enter the commands into a Python interpreter. Use of IPython or Jupyter to run Python will make this a more pleasant experience.

0: Load the GSASIIscriptable module

The first step in writing any Python module is to load the modules that will be needed. Here we will need the os and sys modules from the Python standard library and the GSASIIscriptable from the GSAS-II code. The location where the GSAS-II source code is installed must usually be hard-coded into the script as is done in the example below. Note that a common location for this will be os.path.expanduser("~/g2conda/GSASII/"). Thus the script will begin with:
If a ImportError: No module named GSASIIscriptable error occurs, then the wrong directory location has been supplied for the GSAS-II files.

1: Define some other prerequisite code

To simplify this example, we will define the location where files will be written as workdir (this directory must exist, but it may be empty) and the location where the input files for this exercise may be found as datadir. (Note that this directory must have the the following files: "PBSO4.XRA", "INST_XRY.PRM", "PBSO4.CWN", "inst_d1a.prm" and "PbSO4-Wyckoff.cif". These files can be downloaded from https://subversion.xray.aps.anl.gov/pyGSAS/Tutorials/PythonScript/data/)

We also define here a short function to display the weighted profile R factor for every histogram in a project, HistStats. This will be discussed later.

2: Create a GSAS-II project

The first step in creating a GSASIIscriptable script to to create or access a GSAS-II project, which is done with a call to GSASIIscriptable.G2Project(). This can be done in one of two ways, a call with G2sc.G2Project(filename=file) creates a new (empty) project, while a call with G2sc.G2Project(gpxfile=file) opens and reads an existing project (.gpx) file. Both return a Project wrapper object that is used to access a number of methods and variables. Note that GSASIIscriptable.G2Project() can read .gpx files written by both previous GSASIIscriptable.G2Project() runs or the GSAS-II GUI.

In this example, we create a new project by using the filename=file argument. Note that this will not actually create a file until a save operation is done.

3: Add Histograms and Phase to the GSAS-II project

To add two powder diffraction datasets (histograms) to the project we use the project.add_powder_histogram() method in the GSASIIscriptable class. The two arguments to .add_powder_histogram() are the powder dataset and the instrument parameter file. While neither powder data file uses a standard extension, the importer is able to determine the file format anyway. This will not be true for all file formats.

Then we add a phase to the project using project.add_phase(). This specifies a CIF containing the structural information, a name for the phase and specifies that the two histograms are "added" (linked) to the phase.

4: Change the number of refinement cycles

While this is not noted in the original tutorial, this exercise will not complete properly if more variables are added to the refinement without converging the refinement (or at least coming close to convergence) at each refinement step. This is best accomplished by increasing the number of least-squares cycles. There at present is no method in the project object that allows this parameter to be set, so this must be done by finding appropriate parameter dictionary entry.

At this point, the Python command gpx.data.keys() shows the names of the entries in the GSAS-II tree; 'Controls' is the tree item corresponding to the section where Least Squares cycles are set. Command gpx.data['Controls'].keys() shows that all values are located in an entry labeled 'data' and gpx.data['Controls']['data'].keys() shows the entries in this section. Examination of
gpx.data['Controls']['data']['max cyc']shows the value 3, which is default number of cycles. Thus, the number of cycles is changed with this Python command:

5: Set initial variables and refine

In Step 4 of the original tutorial, the refinement is performed with the variables that are on by default. These variables are the two histogram scale factors and three background parameters for each histogram (8 in total). With GSASIIscriptable, the background refinement flags are not turned on by default, so a dictionary must be created to set these histogram variables. This code: Function HistStats() works by using gpx.histograms() to iterate over all histograms in the project, setting hist to each histogram object. Class member hist.name provides the histogram name and method hist.get_wR() looks up the profile R-factor and prints them. The function also writes the final refinement results into the current project file.
The output from this will be:
 Hessian Levenburg-Marquardt SVD refinement on 8 variables:
initial chi^2 9.6912e+06
 Cycle: 0, Time: 1.88s, Chi**2: 6.7609e+05, Lambda: 0.001,  Delta: 0.93
initial chi^2 6.7609e+05
 Cycle: 1, Time: 1.84s, Chi**2: 6.7602e+05, Lambda: 0.001,  Delta: 0.000104
converged
Found 0 SVD zeros
Read from file: /Users/toby/software/G2/Tutorials/PythonScript/step4.bak0.gpx
Save to file  : /Users/toby/software/G2/Tutorials/PythonScript/step4.gpx
GPX file save successful
 Refinement results are in file: /Users/toby/software/G2/Tutorials/PythonScript/step4.lst
 ***** Refinement successful *****
*** profile Rwp, step4.gpx
	PWDR PBSO4.XRA Bank 1: 40.88
	PWDR PBSO4.CWN Bank 1: 18.65
Note that the Rwp values agree with what is expected from the original tutorial.

Note: that there are several equivalent ways to set the histogram parameters using G2PwdrData.set_refinements(), G2Project.set_refinement() or my_project.do_refinements(), as described here. The gpx.do_refinements([refdict0]) statement above could be replaced with:

gpx.set_refinement({"Background": { "no. coeffs": 3, "refine": True }})
gpx.do_refinements([{}])
or
 
for hist in gpx.histograms():
    hist.set_refinements({"Background": { "no. coeffs": 3, "refine": True }})
gpx.do_refinements([{}])

6: Add unit cell parameters to refinement

In Step 5 of the original tutorial, the refinement is performed again, but the unit cell is refined in the phase.
In this case the gpx.set_refinement(refdict1) statement can be replaced with:
phase0.set_refinements({"Cell": True})
Note that it is also possible to combine the refinement in the current and previous section using
gpx.do_refinements([refdict0,refdict1])
but then the results are not saved in separate project files and it is not possible to see the Rwp values after each refinement.

7: Add hydrostatic strain for just one histogram

In Step 6 of the original tutorial, the refinement is performed again after adding Histogram/Phase (HAP) parameters that allow different lattice constants for the first histogram only.
Here we cannot use gpx.do_refinements([refdict2]) because that would turn on refinement of the Hstrain terms for both histograms (if there were more than one phase, it would be all phases and all histograms), but in the above the gpx.set_refinement(...) statement alternately can be replaced with this:
phase0.set_HAP_refinements({"HStrain": True},histograms=[hist1])

8: Add X-ray Sample broadening terms

The next step in the original tutorial is to treat sample broadening by turning on refinement of the "Mustrain" (microstrain) and "Size" (Scherrer broadening) terms using an isotropic (single-term) model. As described in the documentation for Histogram/Phase parameters, type must always be specfied even as in this case where it is not being changed from the existing default. Again, since these parameters are being set only for one histogram, either phase0.set_HAP_refinements() or gpx.set_refinement() must be used.

9: Add Structural and Sample Displacement Parameters

The step 8 in the original tutorial is to treat sample displacement for each histogram/phase. These parameters are different because of the differing diffractometer geometries. We also add refinement of sample parameters using phase0.set_refinements() to set the "X" (coordinate) and "U" (displacement) refinement flags for all atoms. The original tutorial calls for the diffractometer radius to be changed. This parameter cannot be set from any GSASIIscriptable routines, but following a similar process, as before, the location for this setting can be located in the histogram's 'Sample Parameters' section (hist2.data['Sample Parameters']['Gonio. radius']).
Note that the settings provided in the phase0.set_refinements() statement could have been done using the gpx.do_refinements() call, but not those in the histX.set_refinements() calls, since they must be different for each histogram.

10: Change Data Limits; Vary Gaussian Profile Terms

The final step in the original tutorial is to trim the range of data used in the refinement to exclude data where no reflections occur and where the peaks are cut off at high angle. Also, additional parameters are refined here because the supplied instrumental profile parameters are not very accurate descriptions for the datasets that are used here. It is not possible to refine the Lorentzian x-ray instrumental broadening terms, since this broadening is treated as sample broadening. The Lorentzian neutron broadening is negligible.
Note that the final gpx.do_refinements() statement can be replaced with calls to histX.set_refinements() or gpx.set_refinement(), such as
hist1.set_refinements({'Instrument Parameters': ['U', 'V', 'W']})
hist2.set_refinements({'Instrument Parameters': ['U', 'V', 'W']})
gpx.do_refinements([{}])

Last modified: Thu Oct 12 10:34:10 CDT 2017