source: trunk/GSASIIscriptable.py @ 4021

Last change on this file since 4021 was 4021, checked in by toby, 2 years ago

implement filter for screen messages; start to replace print() with G2Print(); scripting changes for PDF

  • Property svn:eol-style set to native
  • Property svn:keywords set to Date Author Revision URL Id
File size: 231.9 KB
Line 
1#!/usr/bin/env python
2# -*- coding: utf-8 -*-
3########### SVN repository information ###################
4# $Date: 2019-06-12 03:05:07 +0000 (Wed, 12 Jun 2019) $
5# $Author: toby $
6# $Revision: 4021 $
7# $URL: trunk/GSASIIscriptable.py $
8# $Id: GSASIIscriptable.py 4021 2019-06-12 03:05:07Z toby $
9########### SVN repository information ###################
10#
11"""
12*GSASIIscriptable: Scripting Interface*
13=======================================
14
15Routines to use an increasing amount of GSAS-II's capabilities from scripts,
16without use of the graphical user interface (GUI). GSASIIscriptable can create and access
17GSAS-II project (.gpx) files and can directly perform image handling and refinements. 
18The module defines wrapper classes (inheriting from :class:`G2ObjectWrapper`) for a growing number
19of data tree items.
20
21GSASIIscriptable can be used in two ways. It offers a command-line mode
22(see :ref:`CommandlineInterface`) that
23provides access a number of features without writing Python scripts
24via shell/batch commands. The more powerful mode of GSASIIscriptable is
25use is through Python scripts that
26call the module's application interface (API), see API summary that follows or the :ref:`API`
27section.
28
29==================================================
30Application Interface (API) Summary
31==================================================
32This section of the documentation provides an overview to API, with full documentation
33in the :ref:`API` section. The typical API use will be with a Python script, such as this:
34
35.. code-block::  python
36
37    from __future__ import division, print_function
38    import os,sys
39    sys.path.insert(0,'/Users/toby/software/G2/GSASII') # needed to "find" GSAS-II modules
40    import GSASIIscriptable as G2sc
41    datadir = "/Users/Scratch/"
42    gpx = G2sc.G2Project(os.path.join(datadir,'test2.gpx'))
43    gpx.histogram(0).add_back_peak(4.5,30000,5000,0)
44    pardict = {'set': {'Sample Parameters': ['Absorption', 'Contrast', 'DisplaceX'],
45                       'Background': {'type': 'chebyschev', 'refine': True,
46                                      'peaks':[[0,True]]}}}
47    gpx.set_refinement(pardict)
48
49GSASIIscriptable can be used to setup and start sequential refinements. This script
50is used to take the single-dataset fit at the end of Step 1 of the
51`Sequential Refinement tutorial <https://subversion.xray.aps.anl.gov/pyGSAS/Tutorials/SeqRefine/SequentialTutorial.htm>`
52and turn off refinement flags, add histograms and setup the sequential fit, which is then run:
53
54.. code-block::  python
55
56    import os,sys,glob
57    sys.path.insert(0,'/Users/toby/software/G2/GSASII')
58    import GSASIIscriptable as G2sc
59    datadir = os.path.expanduser("~/Scratch/SeqTut2019Mar")
60    PathWrap = lambda fil: os.path.join(datadir,fil)
61    # load and rename project
62    gpx = G2sc.G2Project(PathWrap('7Konly.gpx'))
63    gpx.save(PathWrap('SeqRef.gpx'))
64    # turn off some variables; turn on Dijs
65    for p in gpx.phases():
66        p.set_refinements({"Cell": False})
67    gpx.phase(0).set_HAP_refinements(
68        {'Scale': False,
69         "Size": {'type':'isotropic', 'refine': False},
70         "Mustrain": {'type':'uniaxial', 'refine': False},
71         "HStrain":True,})
72    gpx.phase(1).set_HAP_refinements({'Scale': False})
73    gpx.histogram(0).clear_refinements({'Background':False,
74                     'Sample Parameters':['DisplaceX'],})
75    gpx.histogram(0).ref_back_peak(0,[])
76    gpx.phase(1).set_HAP_refinements({"HStrain":(1,1,1,0)})
77    for fil in sorted(glob.glob(PathWrap('*.fxye'))): # load in remaining fxye files
78        if '00' in fil: continue
79        gpx.add_powder_histogram(fil, PathWrap('OH_00.prm'), fmthint="GSAS powder",phases='all')
80    # copy HAP values, background, instrument params. & limits, not sample params.
81    gpx.copyHistParms(0,'all',['b','i','l'])
82    for p in gpx.phases(): p.copyHAPvalues(0,'all')
83    # setup and launch sequential fit
84    gpx.set_Controls('sequential',gpx.histograms())
85    gpx.set_Controls('cycles',10)
86    gpx.set_Controls('seqCopy',True)
87    gpx.refine() 
88
89Most functionality is provided via the objects and methods described in this section.
90
91---------------------
92:class:`G2Project`
93---------------------
94
95  All GSASIIscriptable scripts will need to create a :class:`G2Project` object
96  either for a new GSAS-II project or to read in an existing project (.gpx) file.
97  The most commonly used routines in this object are:
98
99.. tabularcolumns:: |l|p{3.5in}|
100
101==================================================    ===============================================================================================================
102method                                                Use
103==================================================    ===============================================================================================================
104:meth:`G2Project.save`                                Writes the current project to disk.
105
106:meth:`G2Project.add_powder_histogram`                Used to read in powder diffraction data into a project file.
107
108:meth:`G2Project.add_simulated_powder_histogram`      Defines a "dummy" powder diffraction data that will be simulated after a refinement step.
109
110:meth:`G2Project.add_image`                           Reads in an image into a project.
111
112:meth:`G2Project.add_phase`                           Adds a phase to a project
113
114:meth:`G2Project.add_PDF`                             Adds a PDF entry to a project (does not compute it)
115
116:meth:`G2Project.histograms`                          Provides a list of histograms in the current project, as :class:`G2PwdrData` objects
117
118:meth:`G2Project.phases`                              Provides a list of phases defined in the current project, as :class:`G2Phase` objects
119
120:meth:`G2Project.images`                              Provides a list of images in the current project, as :class:`G2Image` objects
121
122:meth:`G2Project.pdfs`                                Provides a list of PDFs in the current project, as :class:`G2PDF` objects
123
124:meth:`G2Project.seqref`                              Returns a :class:`G2SeqRefRes` object if there are Sequential Refinement results
125
126:meth:`G2Project.do_refinements`                      This is passed a list of dictionaries, where each dict defines a refinement step.
127                                                      Passing a list with a single empty dict initiates a refinement with the current
128                                                      parameters and flags. A refinement dict sets up a single refinement step
129                                                      (as described in :ref:`Project_dicts`). Also see :ref:`Refinement_recipe`.
130
131:meth:`G2Project.set_refinement`                      This is passed a single dict which is used to set parameters and flags.
132                                                      These actions can be performed also in :meth:`G2Project.do_refinements`.
133:meth:`G2Project.get_Variable`                        Retrieves the value and esd for a parameter
134:meth:`G2Project.get_Covariance`                      Retrieves values and covariance for a set of refined parameters
135:meth:`G2Project.set_Controls`                        Set overall GSAS-II control settings such as number of cycles and to set up a sequential
136                                                      fit. (Also see :meth:`G2Project.get_Controls` to read values.)
137==================================================    ===============================================================================================================
138
139---------------------
140:class:`G2Phase`
141---------------------
142
143  Another common object in GSASIIscriptable scripts is :class:`G2Phase`, used to encapsulate each phase in a project, with commonly used methods:
144
145.. tabularcolumns:: |l|p{3.5in}|
146
147==================================================    ===============================================================================================================
148method                                                Use
149==================================================    ===============================================================================================================
150:meth:`G2Phase.set_refinements`                       Provides a mechanism to set values and refinement flags for the phase. See :ref:`Phase_parameters_table`
151                                                      for more details. This information also can be supplied within a call to :meth:`G2Project.do_refinements`
152                                                      or :meth:`G2Project.set_refinement`.
153:meth:`G2Phase.clear_refinements`                     Unsets refinement flags for the phase.
154:meth:`G2Phase.set_HAP_refinements`                   Provides a mechanism to set values and refinement flags for parameters specific to both this phase and
155                                                      one of its histograms. See :ref:`HAP_parameters_table`. This information also can be supplied within
156                                                      a call to :meth:`G2Project.do_refinements` or :meth:`G2Project.set_refinement`.
157:meth:`G2Phase.clear_HAP_refinements`                 Clears refinement flags specific to both this phase and one of its histograms.
158:meth:`G2Phase.getHAPvalues`                          Returns values of parameters specific to both this phase and one of its histograms.
159:meth:`G2Phase.copyHAPvalues`                         Copies HAP settings between from one phase/histogram and to other histograms in same phase.
160:meth:`G2Phase.atoms`                                 Returns a list of atoms in the phase
161:meth:`G2Phase.atom`                                  Returns an atom from its label
162:meth:`G2Phase.histograms`                            Returns a list of histograms linked to the phase
163:meth:`G2Phase.get_cell`                              Returns unit cell parameters (also see :meth:`G2Phase.get_cell_and_esd`)
164:meth:`G2Phase.export_CIF`                            Writes a CIF for the phase
165==================================================    ===============================================================================================================
166
167---------------------
168:class:`G2PwdrData`
169---------------------
170
171  Another common object in GSASIIscriptable scripts is :class:`G2PwdrData`, which encapsulate each powder diffraction histogram in a project, with commonly used methods:
172
173.. tabularcolumns:: |l|p{3.5in}|
174
175==================================================    ===============================================================================================================
176method                                                Use
177==================================================    ===============================================================================================================
178:meth:`G2PwdrData.set_refinements`                    Provides a mechanism to set values and refinement flags for the powder histogram. See
179                                                      :ref:`Histogram_parameters_table` for details. 
180:meth:`G2PwdrData.clear_refinements`                  Unsets refinement flags for the the powder histogram.
181:meth:`G2PwdrData.residuals`                          Reports R-factors etc. for the the powder histogram (also see :meth:`G2PwdrData.get_wR`)
182:meth:`G2PwdrData.add_back_peak`                      Adds a background peak to the histogram. Also see :meth:`G2PwdrData.del_back_peak` and
183                                                      :meth:`G2PwdrData.ref_back_peak`.
184:meth:`G2PwdrData.fit_fixed_points`                   Fits background to the specified fixed points.
185:meth:`G2PwdrData.getdata`                            Provides access to the diffraction data associated with the histogram.
186:meth:`G2PwdrData.Export`                             Writes the diffraction data into a file
187==================================================    ===============================================================================================================
188
189---------------------
190:class:`G2Image`
191---------------------
192
193  When working with images, there will be a :class:`G2Image` object for each image (also see :meth:`G2Project.add_image`  and :meth:`G2Project.images`).
194
195.. tabularcolumns:: |l|p{3.5in}|
196
197==================================================    ===============================================================================================================
198method                                                Use
199==================================================    ===============================================================================================================
200:meth:`G2Image.Recalibrate`                           Invokes a recalibration fit starting from the current Image Controls calibration coefficients.
201:meth:`G2Image.Integrate`                             Invokes an image integration All parameters Image Controls will have previously been set.
202:meth:`G2Image.setControl`                            Set an Image Controls parameter in the current image.
203:meth:`G2Image.getControl`                            Return an Image Controls parameter in the current image.
204:meth:`G2Image.findControl`                           Get the names of Image Controls parameters.
205:meth:`G2Image.loadControls`                          Load controls from a .imctrl file (also see :meth:`G2Image.saveControls`).
206:meth:`G2Image.loadMasks`                             Load masks from a .immask file.
207:meth:`G2Image.setVary`                               Set a refinement flag for Image Controls parameter in the current image. (Also see :meth:`G2Image.getVary`)
208:meth:`G2Image.setCalibrant`                          Set a calibrant type (or show choices) for the current image.
209:meth:`G2Image.setControlFile`                        Set a image to be used as a background/dark/gain map image.
210==================================================    ===============================================================================================================
211
212
213---------------------
214:class:`G2PDF`
215---------------------
216
217  To work with PDF entries, object :class:`G2PDF`, encapsulates a PDF entry with methods:
218
219.. tabularcolumns:: |l|p{3.5in}|
220
221==================================================    ===============================================================================================================
222method                                                Use
223==================================================    ===============================================================================================================
224:meth:`G2PDF.export`                                   Used to write G(r), etc. as a file
225:meth:`G2PDF.calculate`                                Computes the PDF using parameters in the object
226:meth:`G2PDF.optimize`                                 Optimizes selected PDF parameters
227:meth:`G2PDF.set_background`                           Sets the histograms used for sample background, container, etc.
228:meth:`G2PDF.set_formula`                              Sets the chemical formula for the sample
229==================================================    ===============================================================================================================
230
231---------------------
232:class:`G2SeqRefRes`
233---------------------
234
235  To work with Sequential Refinement results, object :class:`G2SeqRefRes`, encapsulates the sequential refinement table with methods:
236
237.. tabularcolumns:: |l|p{3.5in}|
238
239==================================================    ===============================================================================================================
240method                                                Use
241==================================================    ===============================================================================================================
242:meth:`G2SeqRefRes.histograms`                         Provides a list of histograms used in the Sequential Refinement
243:meth:`G2SeqRefRes.get_cell_and_esd`                   Returns cell dimensions and standard uncertainies for a phase and histogram from the Sequential Refinement
244:meth:`G2SeqRefRes.get_Variable`                       Retrieves the value and esd for a parameter from a particular histogram in the Sequential Refinement
245:meth:`G2SeqRefRes.get_Covariance`                     Retrieves values and covariance for a set of refined parameters for a particular histogram
246==================================================    ===============================================================================================================
247
248----------------------
249:class:`G2AtomRecord`
250----------------------
251
252  When working with phases, :class:`G2AtomRecord` objects provide access to the contents of each atom in a phase. This provides access to "properties" that can be
253  used to get values of much of the atoms associated settings: label, type, refinement_flags, coordinates, occupancy, ranId, adp_flag, and uiso. In addition,
254  refinement_flags, occupancy and uiso can be used to set values. See the :class:`G2AtomRecord` docs and source code.
255
256.. _Refinement_dicts:
257
258=====================
259Refinement parameters
260=====================
261While scripts can be written that setup refinements by changing individual parameters
262through calls to the methods associated with objects that wrap each data tree item,
263many of these actions can be combined into fairly complex dict structures to conduct refinement
264steps. Use of these dicts is required with the :ref:`CommandlineInterface`. This section of the
265documentation describes these dicts.
266
267.. _Project_dicts:
268
269-----------------------------
270Project-level Parameter Dict
271-----------------------------
272
273As noted below (:ref:`Refinement_parameters_kinds`), there are three types of refinement parameters,
274which can be accessed individually by the objects that encapsulate individual phases and histograms
275but it will often be simplest to create a composite dictionary
276that is used at the project-level. A dict is created with keys
277"set" and "clear" that can be supplied to :meth:`G2Project.set_refinement`
278(or :meth:`G2Project.do_refinements`, see :ref:`Refinement_recipe` below) that will
279determine parameter values and will determine which parameters will be refined.
280
281The specific keys and subkeys that can be used are defined in tables
282:ref:`Histogram_parameters_table`, :ref:`Phase_parameters_table` and :ref:`HAP_parameters_table`.
283
284Note that optionally a list of histograms and/or phases can be supplied in the call to
285:meth:`G2Project.set_refinement`, but if not specified, the default is to use all defined
286phases and histograms.
287
288As an example:
289
290.. code-block::  python
291
292    pardict = {'set': { 'Limits': [0.8, 12.0],
293                       'Sample Parameters': ['Absorption', 'Contrast', 'DisplaceX'],
294                       'Background': {'type': 'chebyschev', 'refine': True,
295                                      'peaks':[[0,True],[1,1,1]] }},
296              'clear': {'Instrument Parameters': ['U', 'V', 'W']}}
297    my_project.set_refinement(pardict)
298   
299.. _Refinement_recipe:
300   
301------------------------
302Refinement recipe
303------------------------
304Building on the :ref:`Project_dicts`,
305it is possible to specify a sequence of refinement actions as a list of
306these dicts and supplying this list
307as an argument to :meth:`G2Project.do_refinements`.
308
309As an example, this code performs the same actions as in the example in the section above:
310
311.. code-block::  python
312   
313    pardict = {'set': { 'Limits': [0.8, 12.0],
314                       'Sample Parameters': ['Absorption', 'Contrast', 'DisplaceX'],
315                       'Background': {'type': 'chebyschev', 'refine': True}},
316              'clear': {'Instrument Parameters': ['U', 'V', 'W']}}
317    my_project.do_refinements([pardict])
318
319However, in addition to setting a number of parameters, this example will perform a refinement as well,
320after setting the parameters. More than one refinement can be performed by including more
321than one dict in the list.
322
323In this example, two refinement steps will be performed:
324
325.. code-block::  python
326
327    my_project.do_refinements([pardict,pardict1])
328
329
330The keys defined in the following table
331may be used in a dict supplied to :meth:`G2Project.do_refinements`. Note that keys ``histograms``
332and ``phases`` are used to limit actions to specific sets of parameters within the project.
333
334========== ============================================================================
335key         explanation
336========== ============================================================================
337set                    Specifies a dict with keys and subkeys as described in the
338                       :ref:`Refinement_parameters_fmt` section. Items listed here
339                       will be set to be refined.
340clear                  Specifies a dict, as above for set, except that parameters are
341                       cleared and thus will not be refined.
342once                   Specifies a dict as above for set, except that parameters are
343                       set for the next cycle of refinement and are cleared once the
344                       refinement step is completed.
345skip                   Normally, once parameters are processed with a set/clear/once
346                       action(s), a refinement is started. If skip is defined as True
347                       (or any other value) the refinement step is not performed.
348output                 If a file name is specified for output is will be used to save
349                       the current refinement.
350histograms             Should contain a list of histogram(s) to be used for the
351                       set/clear/once action(s) on :ref:`Histogram_parameters_table` or
352                       :ref:`HAP_parameters_table`. Note that this will be
353                       ignored for :ref:`Phase_parameters_table`. Histograms may be
354                       specified as a list of strings [('PWDR ...'),...], indices
355                       [0,1,2] or as list of objects [hist1, hist2].
356phases                 Should contain a list of phase(s) to be used for the
357                       set/clear/once action(s) on :ref:`Phase_parameters_table` or
358                       :ref:`HAP_parameters_table`. Note that this will be
359                       ignored for :ref:`Histogram_parameters_table`.
360                       Phases may be specified as a list of strings
361                       [('Phase name'),...], indices [0,1,2] or as list of objects
362                       [phase0, phase2].
363call                   Specifies a function to call after a refinement is completed.
364                       The value supplied can be the object (typically a function)
365                       that will be called or a string that will evaluate (in the
366                       namespace inside :meth:`G2Project.iter_refinements` where
367                       ``self`` references the project.)
368                       Nothing is called if this is not specified.
369callargs               Provides a list of arguments that will be passed to the function
370                       in call (if any). If call is defined and callargs is not, the
371                       current <tt>G2Project</tt> is passed as a single argument.
372========== ============================================================================
373
374An example that performs a series of refinement steps follows:
375
376.. code-block::  python
377
378    reflist = [
379            {"set": { "Limits": { "low": 0.7 },
380                      "Background": { "no. coeffs": 3,
381                                      "refine": True }}},
382            {"set": { "LeBail": True,
383                      "Cell": True }},
384            {"set": { "Sample Parameters": ["DisplaceX"]}},
385            {"set": { "Instrument Parameters": ["U", "V", "W", "X", "Y"]}},
386            {"set": { "Mustrain": { "type": "uniaxial",
387                                    "refine": "equatorial",
388                                    "direction": [0, 0, 1]}}},
389            {"set": { "Mustrain": { "type": "uniaxial",
390                                    "refine": "axial"}}},
391            {"clear": { "LeBail": True},
392             "set": { "Atoms": { "Mn": "X" }}},
393            {"set": { "Atoms": { "O1": "X", "O2": "X" }}},]
394    my_project.do_refinements(reflist)
395   
396
397In this example, a separate refinement step will be performed for each dict in the list. The keyword
398"skip" can be used to specify a dict that should not include a refinement.
399Note that in the second from last refinement step, parameters are both set and cleared.
400   
401.. _Refinement_parameters_kinds:
402
403----------------------------
404Refinement parameter types
405----------------------------
406
407Note that parameters and refinement flags used in GSAS-II fall into three classes:
408
409    * **Histogram**: There will be a set of these for each dataset loaded into a
410      project file. The parameters available depend on the type of histogram
411      (Bragg-Brentano, Single-Crystal, TOF,...). Typical Histogram parameters
412      include the overall scale factor, background, instrument and sample parameters;
413      see the :ref:`Histogram_parameters_table` table for a list of the histogram
414      parameters where access has been provided.
415     
416    * **Phase**: There will be a set of these for each phase loaded into a
417      project file. While some parameters are found in all types of phases,
418      others are only found in certain types (modulated, magnetic, protein...).
419      Typical phase parameters include unit cell lengths and atomic positions; see the
420      :ref:`Phase_parameters_table` table for a list of the phase     
421      parameters where access has been provided.
422     
423    * **Histogram-and-phase** (HAP): There is a set of these for every histogram
424      that is associated with each phase, so that if there are ``N`` phases and ``M``
425      histograms, there can be ``N*M`` total sets of "HAP" parameters sets (fewer if all
426      histograms are not linked to all phases.) Typical HAP parameters include the
427      phase fractions, sample microstrain and crystallite size broadening terms,
428      hydrostatic strain perturbations of the unit cell and preferred orientation
429      values.
430      See the :ref:`HAP_parameters_table` table for the HAP parameters where access has
431      been provided.
432
433.. _Refinement_parameters_fmt:
434
435=================================
436Specifying Refinement Parameters
437=================================
438
439Refinement parameter values and flags to turn refinement on and off are specified within dictionaries,
440where the details of these dicts are organized depends on the
441type of parameter (see :ref:`Refinement_parameters_kinds`), with a different set
442of keys (as described below) for each of the three types of parameters.
443
444.. _Histogram_parameters_table:
445
446--------------------
447Histogram parameters
448--------------------
449
450This table describes the dictionaries supplied to :func:`G2PwdrData.set_refinements`
451and :func:`G2PwdrData.clear_refinements`. As an example,
452
453.. code-block::  python
454
455   hist.set_refinements({"Background": {"no.coeffs": 3, "refine": True},
456                         "Sample Parameters": ["Scale"],
457                         "Limits": [10000, 40000]})
458
459With :meth:`G2Project.do_refinements`, these parameters should be placed inside a dict with a key
460``set``, ``clear``, or ``once``. Values will be set for all histograms, unless the ``histograms``
461key is used to define specific histograms. As an example:
462
463.. code-block::  python
464
465  gsas_proj.do_refinements([
466      {'set': {
467          'Background': {'no.coeffs': 3, 'refine': True},
468          'Sample Parameters': ['Scale'],
469          'Limits': [10000, 40000]},
470      'histograms': [1,2]}
471                            ])
472
473Note that below in the Instrument Parameters section,
474related profile parameters (such as U and V) are grouped together but
475separated by commas to save space in the table.
476
477.. tabularcolumns:: |l|l|p{3.5in}|
478
479===================== ====================  =================================================
480key                   subkey                explanation
481===================== ====================  =================================================
482Limits                                      The range of 2-theta (degrees) or TOF (in
483                                            microsec) range of values to use. Can
484                                            be either a dictionary of 'low' and/or 'high',
485                                            or a list of 2 items [low, high]
486\                     low                   Sets the low limit
487\                     high                  Sets the high limit
488
489Sample Parameters                           Should be provided as a **list** of subkeys
490                                            to set or clear, e.g. ['DisplaceX', 'Scale']
491\                     Absorption
492\                     Contrast
493\                     DisplaceX             Sample displacement along the X direction
494\                     DisplaceY             Sample displacement along the Y direction
495\                     Scale                 Histogram Scale factor
496
497Background                                  Sample background. Value will be a dict or
498                                            a boolean. If True or False, the refine
499                                            parameter for background is set to that.
500                                            Note that background peaks are not handled
501                                            via this; see
502                                            :meth:`G2PwdrData.ref_back_peak` instead.
503                                            When value is a dict,
504                                            supply any of the following keys:
505\                     type                  The background model, e.g. 'chebyschev'
506\                     refine                The value of the refine flag, boolean
507\                     no. coeffs            Number of coefficients to use, integer
508\                     coeffs                List of floats, literal values for background
509\                     FixedPoints           List of (2-theta, intensity) values for fixed points
510\                     fit fixed points      If True, triggers a fit to the fixed points to
511                                            be calculated. It is calculated when this key is
512                                            detected, regardless of calls to refine.
513                      peaks                 Specifies a set of flags for refining
514                                            background peaks as a nested list. There may
515                                            be an item for each defined background peak
516                                            (or fewer) and each item is a list with the flag
517                                            values for pos,int,sig & gam (fewer than 4 values
518                                            are allowed).
519
520Instrument Parameters                       As in Sample Paramters, provide as a **list** of
521                                            subkeys to
522                                            set or clear, e.g. ['X', 'Y', 'Zero', 'SH/L']
523\                     U, V, W               Gaussian peak profile terms
524\                     X, Y, Z               Lorentzian peak profile terms
525\                     alpha, beta-0,        TOF profile terms
526                      beta-1, beta-q,
527\                     sig-0, sig-1,         TOF profile terms
528                      sig-2, sig-q
529\                     difA, difB, difC      TOF Calibration constants
530\                     Zero                  Zero shift
531\                     SH/L                  Finger-Cox-Jephcoat low-angle peak asymmetry
532\                     Polariz.              Polarization parameter
533\                     Lam                   Lambda, the incident wavelength
534===================== ====================  =================================================
535
536.. _Phase_parameters_table:
537
538----------------
539Phase parameters
540----------------
541
542This table describes the dictionaries supplied to :func:`G2Phase.set_refinements`
543and :func:`G2Phase.clear_refinements`. With :meth:`G2Project.do_refinements`,
544these parameters should be placed inside a dict with a key
545``set``, ``clear``, or ``once``. Values will be set for all phases, unless the ``phases``
546key is used to define specific phase(s).
547
548
549.. tabularcolumns:: |l|p{4.5in}|
550
551======= ==========================================================
552key                   explanation
553======= ==========================================================
554Cell                  Whether or not to refine the unit cell.
555Atoms                 Dictionary of atoms and refinement flags.
556                      Each key should be an atom label, e.g.
557                      'O3', 'Mn5', and each value should be
558                      a string defining what values to refine.
559                      Values can be any combination of 'F'
560                      for fractional occupancy, 'X' for position,
561                      and 'U' for Debye-Waller factor
562LeBail                Enables LeBail intensity extraction.
563======= ==========================================================
564
565
566.. _HAP_parameters_table:
567
568
569Histogram-and-phase parameters
570------------------------------
571
572This table describes the dictionaries supplied to :func:`G2Phase.set_HAP_refinements`
573and :func:`G2Phase.clear_HAP_refinements`. When supplied to
574:meth:`G2Project.do_refinements`, these parameters should be placed inside a dict with a key
575``set``, ``clear``, or ``once``. Values will be set for all histograms used in each phase,
576unless the ``histograms`` and ``phases`` keys are used to define specific phases and histograms.
577
578.. tabularcolumns:: |l|l|p{3.5in}|
579
580=============  ==========  ============================================================
581key             subkey                 explanation
582=============  ==========  ============================================================
583Babinet                                Should be a **list** of the following
584                                       subkeys. If not, assumes both
585                                       BabA and BabU
586\               BabA
587\               BabU
588Extinction                             Boolean, True to refine.
589HStrain                                Boolean or list/tuple, True to refine all
590                                       appropriate D\ :sub:`ij` terms or False
591                                       to not refine any. If a list/tuple, will
592                                       be a set of True & False values for each
593                                       D\ :sub:`ij` term; number of items must
594                                       match number of terms.
595Mustrain
596\               type                   Mustrain model. One of 'isotropic',
597                                       'uniaxial', or 'generalized'. **Should always
598                                       be included when Mustrain is used.**
599\              direction               For uniaxial only. A list of three
600                                       integers,
601                                       the [hkl] direction of the axis.
602\               refine                 Usually boolean, set to True to refine.
603                                       or False to clear.
604                                       For uniaxial model, can specify a value
605                                       of 'axial' or 'equatorial' to set that flag
606                                       to True or a single
607                                       boolean sets both axial and equatorial.
608Size                                   
609\               type                   Size broadening model. One of 'isotropic',
610                                       'uniaxial', or 'ellipsoid'. **Should always
611                                       be specified when Size is used.**
612\              direction               For uniaxial only. A list of three
613                                       integers,
614                                       the [hkl] direction of the axis.
615\               refine                 Boolean, True to refine.
616\               value                  float, size value in microns
617Pref.Ori.                              Boolean, True to refine
618Show                                   Boolean, True to refine
619Use                                    Boolean, True to refine
620Scale                                  Phase fraction; Boolean, True to refine
621=============  ==========  ============================================================
622
623------------------------
624Histogram/Phase objects
625------------------------
626Each phase and powder histogram in a :class:`G2Project` object has an associated
627object. Parameters within each individual object can be turned on and off by calling
628:meth:`G2PwdrData.set_refinements` or :meth:`G2PwdrData.clear_refinements`
629for histogram parameters;
630:meth:`G2Phase.set_refinements` or :meth:`G2Phase.clear_refinements`
631for phase parameters; and :meth:`G2Phase.set_HAP_refinements` or
632:meth:`G2Phase.clear_HAP_refinements`. As an example, if some_histogram is a histogram object (of type :class:`G2PwdrData`), use this to set parameters in that histogram:
633
634.. code-block::  python
635
636    params = { 'Limits': [0.8, 12.0],
637               'Sample Parameters': ['Absorption', 'Contrast', 'DisplaceX'],
638               'Background': {'type': 'chebyschev', 'refine': True}}
639    some_histogram.set_refinements(params)
640
641Likewise to turn refinement flags on, use code such as this:
642
643.. code-block::  python
644
645    params = { 'Instrument Parameters': ['U', 'V', 'W']}
646    some_histogram.set_refinements(params)
647
648and to turn these refinement flags, off use this (Note that the
649``.clear_refinements()`` methods will usually will turn off refinement even
650if a refinement parameter is set in the dict to True.):
651
652.. code-block::  python
653
654    params = { 'Instrument Parameters': ['U', 'V', 'W']}
655    some_histogram.clear_refinements(params)
656
657For phase parameters, use code such as this:
658   
659.. code-block::  python
660
661    params = { 'LeBail': True, 'Cell': True,
662               'Atoms': { 'Mn1': 'X',
663                          'O3': 'XU',
664                          'V4': 'FXU'}}
665    some_histogram.set_refinements(params)
666
667
668and here is an example for HAP parameters:
669
670.. code-block::  python
671
672    params = { 'Babinet': 'BabA',
673               'Extinction': True,
674               'Mustrain': { 'type': 'uniaxial',
675                             'direction': [0, 0, 1],
676                             'refine': True}}
677    some_phase.set_HAP_refinements(params)
678
679Note that the parameters must match the object type and method (phase vs. histogram vs. HAP).
680
681.. _CommandlineInterface:
682
683=======================================
684GSASIIscriptable Command-line Interface
685=======================================
686
687The routines described above are intended to be called from a Python script, but an
688alternate way to access some of the same functionality is to
689invoke the ``GSASIIscriptable.py`` script from
690the command line usually from within a shell script or batch file. This
691will usually be done with a command such as::
692
693       python <path/>GSASIIscriptable.py <subcommand> <file.gpx> <options>
694
695The following subcommands are defined:
696
697        * create, see :func:`create`
698        * add, see :func:`add`
699        * dump, see :func:`dump`
700        * refine, see :func:`refine`
701        * export, :func:`export`
702        * browse, see :func:`IPyBrowse`
703       
704Run::
705
706   python GSASIIscriptable.py --help
707
708to show the available subcommands, and inspect each subcommand with
709`python GSASIIscriptable.py <subcommand> --help` or see the documentation for each of the above routines.
710
711.. _JsonFormat:
712
713-------------------------
714Parameters in JSON files
715-------------------------
716
717The refine command requires two inputs: an existing GSAS-II project (.gpx) file and
718a JSON format file
719(see `Introducing JSON <http://json.org/>`_) that contains a single dict.
720This dict may have two keys:
721
722refinements:
723  This defines the a set of refinement steps in a JSON representation of a
724  :ref:`Refinement_recipe` list.
725
726code:
727  This optionally defines Python code that will be executed after the project is loaded,
728  but before the refinement is started. This can be used to execute Python code to change
729  parameters that are not accessible via a :ref:`Refinement_recipe` dict (note that the
730  project object is accessed with variable ``proj``) or to define code that will be called
731  later (see key ``call`` in the :ref:`Refinement_recipe` section.)
732   
733JSON website: `Introducing JSON <http://json.org/>`_.
734
735.. _API:
736
737============================================================
738API: Complete Documentation
739============================================================
740
741The large number of classes and modules in this module are described below.
742A script will have one or more G2Project objects using :class:`G2Project` and then
743perform actions such as adding a histogram (method :meth:`G2Project.add_powder_histogram`),
744adding a phase (method :meth:`G2Project.add_phase`),
745or setting parameters and performing a refinement
746(method :meth:`G2Project.do_refinements`).
747
748To change settings within histograms, images and phases, one usually needs to use
749methods inside :class:`G2PwdrData`, :class:`G2Image` or :class:`G2Phase`.
750"""
751
752#============================================================================
753# adding a new object type
754# 1) add a new object class (e.g. G2PDF)
755# 2) add the wrapper into G2Project (e.g. _pdfs, pdf, pdfs)
756# 3) add a new method to add the object into a project (G2Project.add_PDF)
757# 4) add to documentation in section :class:`G2Project`
758# 5) add a new documentation section for the new class
759#============================================================================
760
761from __future__ import division, print_function
762import argparse
763import os.path as ospath
764import datetime as dt
765import sys
766import platform
767if '2' in platform.python_version_tuple()[0]:
768    import cPickle
769    strtypes = (str,unicode)
770else:
771    import pickle as cPickle
772    strtypes = (str,bytes)
773import imp
774import copy
775import os
776import random as ran
777
778import numpy.ma as ma
779import scipy.interpolate as si
780import numpy as np
781import scipy as sp
782
783import GSASIIpath
784GSASIIpath.SetBinaryPath(True)  # for now, this is needed before some of these modules can be imported
785import GSASIIobj as G2obj
786import GSASIIpwd as G2pwd
787import GSASIIstrMain as G2strMain
788#import GSASIIIO as G2IO
789import GSASIIstrIO as G2strIO
790import GSASIIspc as G2spc
791import GSASIIElem as G2elem
792import GSASIIfiles as G2fil
793import GSASIIimage as G2img
794
795# Delay imports to not slow down small scripts that don't need them
796Readers = {'Pwdr':[], 'Phase':[], 'Image':[]}
797'''Readers by reader type'''
798exportersByExtension = {}
799'''Specifies the list of extensions that are supported for Powder data export'''
800npsind = lambda x: np.sin(x*np.pi/180.)
801
802def SetPrintLevel(level):
803    '''Set the level of output from calls to :func:`GSASIIfiles.G2Print`,
804    which should be used in place of print() where possible. This is a
805    wrapper for :func:`GSASIIfiles.G2SetPrintLevel` so that this routine is
806    documented here.
807   
808    :param str level: a string used to set the print level, which may be
809    'all', 'warn', 'error' or 'none'.
810    Note that capitalization and extra letters in level are ignored, so
811    'Warn', 'warnings', etc. will all set the mode to 'warn'
812    '''
813    G2fil.G2SetPrintLevel(level)
814    global printLevel
815    for mode in  'all', 'warn', 'error', 'none':
816        if mode in level.lower():
817            printLevel = mode
818            return
819       
820def LoadG2fil():
821    '''Setup GSAS-II importers.
822    Delay importing this module when possible, it is slow.
823    Multiple calls are not. Only the first does anything.
824    '''
825    if len(Readers['Pwdr']) > 0: return
826    # initialize imports
827    Readers['Pwdr'] = G2fil.LoadImportRoutines("pwd", "Powder_Data")
828    Readers['Phase'] = G2fil.LoadImportRoutines("phase", "Phase")
829    Readers['Image'] = G2fil.LoadImportRoutines("img", "Image")
830
831    # initialize exports
832    for obj in G2fil.LoadExportRoutines(None):
833        try:
834            obj.Writer
835        except AttributeError:
836            continue
837        for typ in obj.exporttype:
838            if typ not in exportersByExtension:
839                exportersByExtension[typ] = {obj.extension:obj}
840            else:
841                exportersByExtension[typ][obj.extension] = obj
842
843def LoadDictFromProjFile(ProjFile):
844    '''Read a GSAS-II project file and load items to dictionary
845   
846    :param str ProjFile: GSAS-II project (name.gpx) full file name
847    :returns: Project,nameList, where
848
849      * Project (dict) is a representation of gpx file following the GSAS-II tree structure
850        for each item: key = tree name (e.g. 'Controls','Restraints',etc.), data is dict
851        data dict = {'data':item data whch may be list, dict or None,'subitems':subdata (if any)}
852      * nameList (list) has names of main tree entries & subentries used to reconstruct project file
853
854    Example for fap.gpx::
855
856      Project = {                 #NB:dict order is not tree order
857        'Phases':{'data':None,'fap':{phase dict}},
858        'PWDR FAP.XRA Bank 1':{'data':[histogram data list],'Comments':comments,'Limits':limits, etc},
859        'Rigid bodies':{'data': {rigid body dict}},
860        'Covariance':{'data':{covariance data dict}},
861        'Controls':{'data':{controls data dict}},
862        'Notebook':{'data':[notebook list]},
863        'Restraints':{'data':{restraint data dict}},
864        'Constraints':{'data':{constraint data dict}}]
865        }
866      nameList = [                #NB: reproduces tree order
867        ['Notebook',],
868        ['Controls',],
869        ['Covariance',],
870        ['Constraints',],
871        ['Restraints',],
872        ['Rigid bodies',],
873        ['PWDR FAP.XRA Bank 1',
874             'Comments',
875             'Limits',
876             'Background',
877             'Instrument Parameters',
878             'Sample Parameters',
879             'Peak List',
880             'Index Peak List',
881             'Unit Cells List',
882             'Reflection Lists'],
883        ['Phases', 'fap']
884        ]
885    '''
886    # Let IOError be thrown if file does not exist
887    if not ospath.exists(ProjFile):
888        G2fil.G2Print ('\n*** Error attempt to open project file that does not exist: \n    {}'.
889                   format(ProjFile))
890        raise IOError('GPX file {} does not exist'.format(ProjFile))
891    try:
892        Project, nameList = G2strIO.GetFullGPX(ProjFile)
893    except Exception as msg:
894        raise IOError(msg)
895    return Project,nameList
896
897def SaveDictToProjFile(Project,nameList,ProjFile):
898    '''Save a GSAS-II project file from dictionary/nameList created by LoadDictFromProjFile
899
900    :param dict Project: representation of gpx file following the GSAS-II
901        tree structure as described for LoadDictFromProjFile
902    :param list nameList: names of main tree entries & subentries used to reconstruct project file
903    :param str ProjFile: full file name for output project.gpx file (including extension)
904    '''
905    file = open(ProjFile,'wb')
906    try:
907        for name in nameList:
908            data = []
909            item = Project[name[0]]
910            data.append([name[0],item['data']])
911            for item2 in name[1:]:
912                data.append([item2,item[item2]])
913            cPickle.dump(data,file,1)
914    finally:
915        file.close()
916    G2fil.G2Print('gpx file saved as %s'%ProjFile)
917
918# def ImportPowder(reader,filename):
919#     '''Use a reader to import a powder diffraction data file
920
921#     :param str reader: a scriptable reader
922#     :param str filename: full name of powder data file; can be "multi-Bank" data
923
924#     :returns: list rdlist: list of reader objects containing powder data, one for each
925#         "Bank" of data encountered in file. Items in reader object of interest are:
926
927#           * rd.comments: list of str: comments found on powder file
928#           * rd.dnames: list of str: data nammes suitable for use in GSASII data tree NB: duplicated in all rd entries in rdlist
929#           * rd.powderdata: list of numpy arrays: pos,int,wt,zeros,zeros,zeros as needed for a PWDR entry in  GSASII data tree.
930#     '''
931#     rdfile,rdpath,descr = imp.find_module(reader)
932#     rdclass = imp.load_module(reader,rdfile,rdpath,descr)
933#     rd = rdclass.GSAS_ReaderClass()
934#     if not rd.scriptable:
935#         G2fil.G2Print(u'**** ERROR: '+reader+u' is not a scriptable reader')
936#         return None
937#     rdlist = []
938#     if rd.ContentsValidator(filename):
939#         repeat = True
940#         rdbuffer = {} # create temporary storage for file reader
941#         block = 0
942#         while repeat: # loop if the reader asks for another pass on the file
943#             block += 1
944#             repeat = False
945#             rd.objname = ospath.basename(filename)
946#             flag = rd.Reader(filename,None,buffer=rdbuffer,blocknum=block,)
947#             if flag:
948#                 rdlist.append(copy.deepcopy(rd)) # save the result before it is written over
949#                 if rd.repeat:
950#                     repeat = True
951#         return rdlist
952#     G2fil.G2Print(rd.errors)
953#     return None
954
955def SetDefaultDData(dType,histoName,NShkl=0,NDij=0):
956    '''Create an initial Histogram dictionary
957
958    Author: Jackson O'Donnell (jacksonhodonnell .at. gmail.com)
959    '''
960    if dType in ['SXC','SNC']:
961        return {'Histogram':histoName,'Show':False,'Scale':[1.0,True],
962            'Babinet':{'BabA':[0.0,False],'BabU':[0.0,False]},
963            'Extinction':['Lorentzian','None', {'Tbar':0.1,'Cos2TM':0.955,
964            'Eg':[1.e-10,False],'Es':[1.e-10,False],'Ep':[1.e-10,False]}],
965            'Flack':[0.0,False]}
966    elif dType == 'SNT':
967        return {'Histogram':histoName,'Show':False,'Scale':[1.0,True],
968            'Babinet':{'BabA':[0.0,False],'BabU':[0.0,False]},
969            'Extinction':['Lorentzian','None', {
970            'Eg':[1.e-10,False],'Es':[1.e-10,False],'Ep':[1.e-10,False]}]}
971    elif 'P' in dType:
972        return {'Histogram':histoName,'Show':False,'Scale':[1.0,False],
973            'Pref.Ori.':['MD',1.0,False,[0,0,1],0,{},[],0.1],
974            'Size':['isotropic',[1.,1.,1.],[False,False,False],[0,0,1],
975                [1.,1.,1.,0.,0.,0.],6*[False,]],
976            'Mustrain':['isotropic',[1000.0,1000.0,1.0],[False,False,False],[0,0,1],
977                NShkl*[0.01,],NShkl*[False,]],
978            'HStrain':[NDij*[0.0,],NDij*[False,]],
979            'Extinction':[0.0,False],'Babinet':{'BabA':[0.0,False],'BabU':[0.0,False]}}
980
981
982def PreSetup(data):
983    '''Create part of an initial (empty) phase dictionary
984
985    from GSASIIphsGUI.py, near end of UpdatePhaseData
986
987    Author: Jackson O'Donnell (jacksonhodonnell .at. gmail.com)
988    '''
989    if 'RBModels' not in data:
990        data['RBModels'] = {}
991    if 'MCSA' not in data:
992        data['MCSA'] = {'Models':[{'Type':'MD','Coef':[1.0,False,[.8,1.2],],'axis':[0,0,1]}],'Results':[],'AtInfo':{}}
993    if 'dict' in str(type(data['MCSA']['Results'])):
994        data['MCSA']['Results'] = []
995    if 'Modulated' not in data['General']:
996        data['General']['Modulated'] = False
997#    if 'modulated' in data['General']['Type']:
998#        data['General']['Modulated'] = True
999#        data['General']['Type'] = 'nuclear'
1000
1001
1002def SetupGeneral(data, dirname):
1003    """Helps initialize phase data.
1004
1005    From GSASIIphsGui.py, function of the same name. Minor changes for imports etc.
1006
1007    Author: Jackson O'Donnell (jacksonhodonnell .at. gmail.com)
1008    """
1009    mapDefault = {'MapType':'','RefList':'','Resolution':0.5,'Show bonds':True,
1010                'rho':[],'rhoMax':0.,'mapSize':10.0,'cutOff':50.,'Flip':False}
1011    generalData = data['General']
1012    atomData = data['Atoms']
1013    generalData['AtomTypes'] = []
1014    generalData['Isotopes'] = {}
1015
1016    if 'Isotope' not in generalData:
1017        generalData['Isotope'] = {}
1018    if 'Data plot type' not in generalData:
1019        generalData['Data plot type'] = 'Mustrain'
1020    if 'POhkl' not in generalData:
1021        generalData['POhkl'] = [0,0,1]
1022    if 'Map' not in generalData:
1023        generalData['Map'] = mapDefault.copy()
1024    if 'Flip' not in generalData:
1025        generalData['Flip'] = {'RefList':'','Resolution':0.5,'Norm element':'None',
1026            'k-factor':0.1,'k-Max':20.,}
1027    if 'testHKL' not in generalData['Flip']:
1028        generalData['Flip']['testHKL'] = [[0,0,2],[2,0,0],[1,1,1],[0,2,0],[1,2,3]]
1029    if 'doPawley' not in generalData:
1030        generalData['doPawley'] = False     #ToDo: change to ''
1031    if 'Pawley dmin' not in generalData:
1032        generalData['Pawley dmin'] = 1.0
1033    if 'Pawley neg wt' not in generalData:
1034        generalData['Pawley neg wt'] = 0.0
1035    if 'Algolrithm' in generalData.get('MCSA controls',{}) or \
1036        'MCSA controls' not in generalData:
1037        generalData['MCSA controls'] = {'Data source':'','Annealing':[50.,0.001,50],
1038        'dmin':2.0,'Algorithm':'log','Jump coeff':[0.95,0.5],'boltzmann':1.0,
1039        'fast parms':[1.0,1.0,1.0],'log slope':0.9,'Cycles':1,'Results':[],'newDmin':True}
1040    if 'AtomPtrs' not in generalData:
1041        generalData['AtomPtrs'] = [3,1,7,9]
1042        if generalData['Type'] == 'macromolecular':
1043            generalData['AtomPtrs'] = [6,4,10,12]
1044        elif generalData['Type'] == 'magnetic':
1045            generalData['AtomPtrs'] = [3,1,10,12]
1046    if generalData['Modulated']:
1047        generalData['Type'] = 'nuclear'
1048        if 'Super' not in generalData:
1049            generalData['Super'] = 1
1050            generalData['SuperVec'] = [[0,0,.1],False,4]
1051            generalData['SSGData'] = {}
1052        if '4DmapData' not in generalData:
1053            generalData['4DmapData'] = mapDefault.copy()
1054            generalData['4DmapData'].update({'MapType':'Fobs'})
1055    if 'Modulated' not in generalData:
1056        generalData['Modulated'] = False
1057    if 'HydIds' not in generalData:
1058        generalData['HydIds'] = {}
1059    cx,ct,cs,cia = generalData['AtomPtrs']
1060    generalData['NoAtoms'] = {}
1061    generalData['BondRadii'] = []
1062    generalData['AngleRadii'] = []
1063    generalData['vdWRadii'] = []
1064    generalData['AtomMass'] = []
1065    generalData['Color'] = []
1066    if generalData['Type'] == 'magnetic':
1067        generalData['MagDmin'] = generalData.get('MagDmin',1.0)
1068        landeg = generalData.get('Lande g',[])
1069    generalData['Mydir'] = dirname
1070    badList = {}
1071    for iat,atom in enumerate(atomData):
1072        atom[ct] = atom[ct].lower().capitalize()              #force to standard form
1073        if generalData['AtomTypes'].count(atom[ct]):
1074            generalData['NoAtoms'][atom[ct]] += atom[cx+3]*float(atom[cs+1])
1075        elif atom[ct] != 'UNK':
1076            Info = G2elem.GetAtomInfo(atom[ct])
1077            if not Info:
1078                if atom[ct] not in badList:
1079                    badList[atom[ct]] = 0
1080                badList[atom[ct]] += 1
1081                atom[ct] = 'UNK'
1082                continue
1083            atom[ct] = Info['Symbol'] # N.B. symbol might be changed by GetAtomInfo
1084            generalData['AtomTypes'].append(atom[ct])
1085            generalData['Z'] = Info['Z']
1086            generalData['Isotopes'][atom[ct]] = Info['Isotopes']
1087            generalData['BondRadii'].append(Info['Drad'])
1088            generalData['AngleRadii'].append(Info['Arad'])
1089            generalData['vdWRadii'].append(Info['Vdrad'])
1090            if atom[ct] in generalData['Isotope']:
1091                if generalData['Isotope'][atom[ct]] not in generalData['Isotopes'][atom[ct]]:
1092                    isotope = list(generalData['Isotopes'][atom[ct]].keys())[-1]
1093                    generalData['Isotope'][atom[ct]] = isotope
1094                generalData['AtomMass'].append(Info['Isotopes'][generalData['Isotope'][atom[ct]]]['Mass'])
1095            else:
1096                generalData['Isotope'][atom[ct]] = 'Nat. Abund.'
1097                if 'Nat. Abund.' not in generalData['Isotopes'][atom[ct]]:
1098                    isotope = list(generalData['Isotopes'][atom[ct]].keys())[-1]
1099                    generalData['Isotope'][atom[ct]] = isotope
1100                generalData['AtomMass'].append(Info['Mass'])
1101            generalData['NoAtoms'][atom[ct]] = atom[cx+3]*float(atom[cs+1])
1102            generalData['Color'].append(Info['Color'])
1103            if generalData['Type'] == 'magnetic':
1104                if len(landeg) < len(generalData['AtomTypes']):
1105                    landeg.append(2.0)
1106    if generalData['Type'] == 'magnetic':
1107        generalData['Lande g'] = landeg[:len(generalData['AtomTypes'])]
1108
1109    if badList:
1110        msg = 'Warning: element symbol(s) not found:'
1111        for key in badList:
1112            msg += '\n\t' + key
1113            if badList[key] > 1:
1114                msg += ' (' + str(badList[key]) + ' times)'
1115        raise G2ScriptException("Phase error:\n" + msg)
1116        # wx.MessageBox(msg,caption='Element symbol error')
1117    F000X = 0.
1118    F000N = 0.
1119    for i,elem in enumerate(generalData['AtomTypes']):
1120        F000X += generalData['NoAtoms'][elem]*generalData['Z']
1121        isotope = generalData['Isotope'][elem]
1122        F000N += generalData['NoAtoms'][elem]*generalData['Isotopes'][elem][isotope]['SL'][0]
1123    generalData['F000X'] = F000X
1124    generalData['F000N'] = F000N
1125    import GSASIImath as G2mth
1126    generalData['Mass'] = G2mth.getMass(generalData)
1127
1128
1129def make_empty_project(author=None, filename=None):
1130    """Creates an dictionary in the style of GSASIIscriptable, for an empty
1131    project.
1132
1133    If no author name or filename is supplied, 'no name' and
1134    <current dir>/test_output.gpx are used , respectively.
1135
1136    Returns: project dictionary, name list
1137
1138    Author: Jackson O'Donnell (jacksonhodonnell .at. gmail.com)
1139    """
1140    if not filename:
1141        filename = 'test_output.gpx'
1142    filename = os.path.abspath(filename)
1143    gsasii_version = str(GSASIIpath.GetVersionNumber())
1144    LoadG2fil()
1145    try:
1146        import matplotlib as mpl
1147        python_library_versions = G2fil.get_python_versions([mpl, np, sp])
1148    except ImportError:
1149        python_library_versions = G2fil.get_python_versions([np, sp])
1150
1151    controls_data = dict(G2obj.DefaultControls)
1152    controls_data['LastSavedAs'] = filename
1153    controls_data['LastSavedUsing'] = gsasii_version
1154    controls_data['PythonVersions'] = python_library_versions
1155    if author:
1156        controls_data['Author'] = author
1157
1158    output = {'Constraints': {'data': {'HAP': [], 'Hist': [], 'Phase': [],
1159                                       'Global': []}},
1160              'Controls': {'data': controls_data},
1161              u'Covariance': {'data': {}},
1162              u'Notebook': {'data': ['']},
1163              u'Restraints': {'data': {}},
1164              u'Rigid bodies': {'data': {'RBIds': {'Residue': [], 'Vector': []},
1165                                'Residue': {'AtInfo': {}},
1166                                'Vector':  {'AtInfo': {}}}}}
1167
1168    names = [[u'Notebook'], [u'Controls'], [u'Covariance'],
1169             [u'Constraints'], [u'Restraints'], [u'Rigid bodies']]
1170
1171    return output, names
1172
1173
1174class G2ImportException(Exception):
1175    pass
1176
1177class G2ScriptException(Exception):
1178    pass
1179
1180def import_generic(filename, readerlist, fmthint=None, bank=None):
1181    """Attempt to import a filename, using a list of reader objects.
1182
1183    Returns the first reader object which worked."""
1184    # Translated from OnImportGeneric method in GSASII.py
1185    primaryReaders, secondaryReaders = [], []
1186    for reader in readerlist:
1187        if fmthint is not None and fmthint not in reader.formatName: continue
1188        flag = reader.ExtensionValidator(filename)
1189        if flag is None:
1190            secondaryReaders.append(reader)
1191        elif flag:
1192            primaryReaders.append(reader)
1193    if not secondaryReaders and not primaryReaders:
1194        raise G2ImportException("Could not read file: ", filename)
1195
1196    with open(filename, 'Ur') as fp:
1197        rd_list = []
1198
1199        for rd in primaryReaders + secondaryReaders:
1200            # Initialize reader
1201            rd.selections = []
1202            if bank is None:
1203                rd.selections = []
1204            else:
1205                rd.selections = [bank-1]
1206            rd.dnames = []
1207            rd.ReInitialize()
1208            # Rewind file
1209            rd.errors = ""
1210            if not rd.ContentsValidator(filename):
1211                # Report error
1212                G2fil.G2Print("Warning: File {} has a validation error, continuing".format(filename))
1213            if len(rd.selections) > 1:
1214                raise G2ImportException("File {} has {} banks. Specify which bank to read with databank param."
1215                                .format(filename,len(rd.selections)))
1216
1217            block = 0
1218            rdbuffer = {}
1219            repeat = True
1220            while repeat:
1221                repeat = False
1222                block += 1
1223                rd.objname = os.path.basename(filename)
1224                try:
1225                    flag = rd.Reader(filename,buffer=rdbuffer, blocknum=block)
1226                except:
1227                    flag = False
1228                if flag:
1229                    # Omitting image loading special cases
1230                    rd.readfilename = filename
1231                    rd_list.append(copy.deepcopy(rd))
1232                    repeat = rd.repeat
1233                else:
1234                    G2fil.G2Print("Warning: {} Reader failed to read {}".format(rd.formatName,filename))
1235            if rd_list:
1236                if rd.warnings:
1237                    G2fil.G2Print("Read warning by", rd.formatName, "reader:",
1238                          rd.warnings)
1239                elif bank is None:
1240                    G2fil.G2Print("{} read by Reader {}"
1241                              .format(filename,rd.formatName))
1242                else:
1243                    G2fil.G2Print("{} block # {} read by Reader {}"
1244                              .format(filename,bank,rd.formatName))
1245                return rd_list
1246    raise G2ImportException("No reader could read file: " + filename)
1247
1248
1249def load_iprms(instfile, reader, bank=None):
1250    """Loads instrument parameters from a file, and edits the
1251    given reader.
1252
1253    Returns a 2-tuple of (Iparm1, Iparm2) parameters
1254    """
1255    LoadG2fil()
1256    ext = os.path.splitext(instfile)[1]
1257
1258    if ext.lower() == '.instprm':
1259        # New GSAS File, load appropriate bank
1260        with open(instfile) as f:
1261            lines = f.readlines()
1262        if bank is None: 
1263            bank = reader.powderentry[2] 
1264        numbanks = reader.numbanks
1265        iparms = G2fil.ReadPowderInstprm(lines, bank, numbanks, reader)
1266        reader.instfile = instfile
1267        reader.instmsg = '{} (G2 fmt) bank {}'.format(instfile,bank)
1268        return iparms
1269    elif ext.lower() not in ('.prm', '.inst', '.ins'):
1270        raise ValueError('Expected .prm file, found: ', instfile)
1271
1272    # It's an old GSAS file, load appropriately
1273    Iparm = {}
1274    with open(instfile, 'Ur') as fp:
1275        for line in fp:
1276            if '#' in line:
1277                continue
1278            Iparm[line[:12]] = line[12:-1]
1279    ibanks = int(Iparm.get('INS   BANK  ', '1').strip())
1280    if bank is not None:
1281        # pull out requested bank # bank from the data, and change the bank to 1
1282        Iparm,IparmC = {},Iparm
1283        for key in IparmC:
1284            if 'INS' not in key[:3]: continue   #skip around rubbish lines in some old iparm
1285            if key[4:6] == "  ":
1286                Iparm[key] = IparmC[key]
1287            elif int(key[4:6].strip()) == bank:
1288                Iparm[key[:4]+' 1'+key[6:]] = IparmC[key]           
1289        reader.instbank = bank
1290    elif ibanks == 1:
1291        reader.instbank = 1
1292    else: 
1293        raise G2ImportException("Instrument parameter file has {} banks, select one with instbank param."
1294                                    .format(ibanks))
1295    reader.powderentry[2] = 1
1296    reader.instfile = instfile
1297    reader.instmsg = '{} bank {}'.format(instfile,reader.instbank)
1298    return G2fil.SetPowderInstParms(Iparm, reader)
1299
1300def load_pwd_from_reader(reader, instprm, existingnames=[],bank=None):
1301    """Loads powder data from a reader object, and assembles it into a GSASII data tree.
1302
1303    :returns: (name, tree) - 2-tuple of the histogram name (str), and data
1304
1305    Author: Jackson O'Donnell (jacksonhodonnell .at. gmail.com)
1306    """
1307    HistName = 'PWDR ' + G2obj.StripUnicode(reader.idstring, '_')
1308    HistName = G2obj.MakeUniqueLabel(HistName, existingnames)
1309
1310    try:
1311        Iparm1, Iparm2 = instprm
1312    except ValueError:
1313        Iparm1, Iparm2 = load_iprms(instprm, reader, bank=bank)
1314        G2fil.G2Print('Instrument parameters read:',reader.instmsg)
1315    Ymin = np.min(reader.powderdata[1])
1316    Ymax = np.max(reader.powderdata[1])
1317    valuesdict = {'wtFactor': 1.0,
1318                  'Dummy': False,
1319                  'ranId': ran.randint(0, sys.maxsize),
1320                  'Offset': [0.0, 0.0], 'delOffset': 0.02*Ymax,
1321                  'refOffset': -0.1*Ymax, 'refDelt': 0.1*Ymax,
1322                  'Yminmax': [Ymin, Ymax]}
1323    reader.Sample['ranId'] = valuesdict['ranId']
1324    if 'T' in Iparm1['Type'][0]:
1325        if not reader.clockWd and reader.GSAS:
1326            reader.powderdata[0] *= 100.        #put back the CW centideg correction
1327
1328    # Ending keys:
1329    # [u'Reflection Lists',
1330    #  u'Limits',
1331    #  'data',
1332    #  u'Index Peak List',
1333    #  u'Comments',
1334    #  u'Unit Cells List',
1335    #  u'Sample Parameters',
1336    #  u'Peak List',
1337    #  u'Background',
1338    #  u'Instrument Parameters']
1339    Tmin = np.min(reader.powderdata[0])
1340    Tmax = np.max(reader.powderdata[0])
1341
1342    default_background = [['chebyschev', False, 3, 1.0, 0.0, 0.0],
1343                          {'nDebye': 0, 'debyeTerms': [], 'nPeaks': 0, 'peaksList': []}]
1344
1345    output_dict = {u'Reflection Lists': {},
1346                   u'Limits': reader.pwdparms.get('Limits', [(Tmin, Tmax), [Tmin, Tmax]]),
1347                   u'data': [valuesdict, reader.powderdata, HistName],
1348                   u'Index Peak List': [[], []],
1349                   u'Comments': reader.comments,
1350                   u'Unit Cells List': [],
1351                   u'Sample Parameters': reader.Sample,
1352                   u'Peak List': {'peaks': [], 'sigDict': {}},
1353                   u'Background': reader.pwdparms.get('Background', default_background),
1354                   u'Instrument Parameters': [Iparm1, Iparm2],
1355                   }
1356
1357    names = [u'Comments',
1358             u'Limits',
1359             u'Background',
1360             u'Instrument Parameters',
1361             u'Sample Parameters',
1362             u'Peak List',
1363             u'Index Peak List',
1364             u'Unit Cells List',
1365             u'Reflection Lists']
1366
1367    # TODO controls?? GSASII.py:1664-7
1368
1369    return HistName, [HistName] + names, output_dict
1370
1371def _deep_copy_into(from_, into):
1372    """Helper function for reloading .gpx file. See G2Project.reload()
1373
1374    :author: Jackson O'Donnell (jacksonhodonnell .at. gmail.com)
1375    """
1376    if isinstance(from_, dict) and isinstance(into, dict):
1377        combined_keys = set(from_.keys()).union(into.keys())
1378        for key in combined_keys:
1379            if key in from_ and key in into:
1380                both_dicts = (isinstance(from_[key], dict)
1381                              and isinstance(into[key], dict))
1382                both_lists = (isinstance(from_[key], list)
1383                              and isinstance(into[key], list))
1384                if both_dicts or both_lists:
1385                    _deep_copy_into(from_[key], into[key])
1386                else:
1387                    into[key] = from_[key]
1388            elif key in from_:
1389                into[key] = from_[key]
1390            else:  # key in into
1391                del into[key]
1392    elif isinstance(from_, list) and isinstance(into, list):
1393        if len(from_) == len(into):
1394            for i in range(len(from_)):
1395                both_dicts = (isinstance(from_[i], dict)
1396                              and isinstance(into[i], dict))
1397                both_lists = (isinstance(from_[i], list)
1398                              and isinstance(into[i], list))
1399                if both_dicts or both_lists:
1400                    _deep_copy_into(from_[i], into[i])
1401                else:
1402                    into[i] = from_[i]
1403        else:
1404            into[:] = from_
1405
1406def _getCorrImage(ImageReaderlist,proj,imageRef):
1407    '''Gets image & applies dark, background & flat background corrections.
1408    based on :func:`GSASIIimgGUI.GetImageZ`. Expected to be for internal
1409    use only.
1410
1411    :param list ImageReaderlist: list of Reader objects for images
1412    :param object ImageReaderlist: list of Reader objects for images
1413    :param imageRef: A reference to the desired image. Either the Image
1414      tree name (str), the image's index (int) or
1415      a image object (:class:`G2Image`)
1416
1417    :return: array sumImg: corrected image for background/dark/flat back
1418    '''
1419    ImgObj = proj.image(imageRef)
1420    Controls = ImgObj.data['Image Controls']
1421    formatName = Controls.get('formatName','')
1422    imagefile = ImgObj.data['data'][1]
1423    ImageTag = None # fix this for multiimage files
1424    sumImg = G2fil.RereadImageData(ImageReaderlist,imagefile,ImageTag=ImageTag,FormatName=formatName)
1425    if sumImg is None:
1426        return []
1427    darkImg = False
1428    if 'dark image' in Controls:
1429        darkImg,darkScale = Controls['dark image']
1430        if darkImg:
1431            dImgObj = proj.image(darkImg)
1432            formatName = dImgObj.data['Image Controls'].get('formatName','')
1433            imagefile = dImgObj.data['data'][1]
1434            ImageTag = None # fix this for multiimage files
1435            darkImage = G2fil.RereadImageData(ImageReaderlist,imagefile,ImageTag=ImageTag,FormatName=formatName)
1436            if darkImg is None:
1437                raise Exception('Error reading dark image {}'.format(imagefile))
1438            sumImg += np.array(darkImage*darkScale,dtype='int32')
1439    if 'background image' in Controls:
1440        backImg,backScale = Controls['background image']           
1441        if backImg:     #ignores any transmission effect in the background image
1442            bImgObj = proj.image(backImg)
1443            formatName = bImgObj.data['Image Controls'].get('formatName','')
1444            imagefile = bImgObj.data['data'][1]
1445            ImageTag = None # fix this for multiimage files
1446            backImage = G2fil.RereadImageData(ImageReaderlist,imagefile,ImageTag=ImageTag,FormatName=formatName)
1447            if backImage is None:
1448                raise Exception('Error reading background image {}'.format(imagefile))
1449            if darkImg:
1450                backImage += np.array(darkImage*darkScale/backScale,dtype='int32')
1451            else:
1452                sumImg += np.array(backImage*backScale,dtype='int32')
1453    if 'Gain map' in Controls:
1454        gainMap = Controls['Gain map']
1455        if gainMap:
1456            gImgObj = proj.image(gainMap)
1457            formatName = gImgObj.data['Image Controls'].get('formatName','')
1458            imagefile = gImgObj.data['data'][1]
1459            ImageTag = None # fix this for multiimage files
1460            GMimage = G2fil.RereadImageData(ImageReaderlist,imagefile,ImageTag=ImageTag,FormatName=formatName)
1461            if GMimage is None:
1462                raise Exception('Error reading Gain map image {}'.format(imagefile))
1463            sumImg = sumImg*GMimage/1000
1464    sumImg -= int(Controls.get('Flat Bkg',0))
1465    Imax = np.max(sumImg)
1466    Controls['range'] = [(0,Imax),[0,Imax]]
1467    return np.asarray(sumImg,dtype='int32')
1468
1469class G2ObjectWrapper(object):
1470    """Base class for all GSAS-II object wrappers.
1471
1472    The underlying GSAS-II format can be accessed as `wrapper.data`. A number
1473    of overrides are implemented so that the wrapper behaves like a dictionary.
1474
1475    Author: Jackson O'Donnell (jacksonhodonnell .at. gmail.com)
1476    """
1477    def __init__(self, datadict):
1478        self.data = datadict
1479
1480    def __getitem__(self, key):
1481        return self.data[key]
1482
1483    def __setitem__(self, key, value):
1484        self.data[key] = value
1485
1486    def __contains__(self, key):
1487        return key in self.data
1488
1489    def get(self, k, d=None):
1490        return self.data.get(k, d)
1491
1492    def keys(self):
1493        return self.data.keys()
1494
1495    def values(self):
1496        return self.data.values()
1497
1498    def items(self):
1499        return self.data.items()
1500
1501
1502class G2Project(G2ObjectWrapper):   
1503    """Represents an entire GSAS-II project.
1504
1505    :param str gpxfile: Existing .gpx file to be loaded. If nonexistent,
1506            creates an empty project.
1507    :param str author: Author's name (not yet implemented)
1508    :param str newgpx: The filename the project should be saved to in
1509            the future. If both newgpx and gpxfile are present, the project is
1510            loaded from the gpxfile, then when saved will be written to newgpx.
1511    :param str filename: Name to be used to save the project. Has same function as
1512            parameter newgpx (do not use both gpxfile and filename). Use of newgpx
1513            is preferred over filename.
1514
1515    There are two ways to initialize this object:
1516
1517    >>> # Load an existing project file
1518    >>> proj = G2Project('filename.gpx')
1519   
1520    >>> # Create a new project
1521    >>> proj = G2Project(newgpx='new_file.gpx')
1522   
1523    Histograms can be accessed easily.
1524
1525    >>> # By name
1526    >>> hist = proj.histogram('PWDR my-histogram-name')
1527   
1528    >>> # Or by index
1529    >>> hist = proj.histogram(0)
1530    >>> assert hist.id == 0
1531   
1532    >>> # Or by random id
1533    >>> assert hist == proj.histogram(hist.ranId)
1534
1535    Phases can be accessed the same way.
1536
1537    >>> phase = proj.phase('name of phase')
1538
1539    New data can also be loaded via :meth:`~G2Project.add_phase` and
1540    :meth:`~G2Project.add_powder_histogram`.
1541
1542    >>> hist = proj.add_powder_histogram('some_data_file.chi',
1543                                         'instrument_parameters.prm')
1544    >>> phase = proj.add_phase('my_phase.cif', histograms=[hist])
1545
1546    Parameters for Rietveld refinement can be turned on and off as well.
1547    See :meth:`~G2Project.set_refinement`, :meth:`~G2Project.clear_refinements`,
1548    :meth:`~G2Project.iter_refinements`, :meth:`~G2Project.do_refinements`.
1549    """
1550    def __init__(self, gpxfile=None, author=None, filename=None, newgpx=None):
1551        if filename is not None and newgpx is not None:
1552            raise G2ScriptException('Do not use filename and newgpx together')
1553        elif newgpx:
1554            filename = newgpx
1555        if not gpxfile:
1556            filename = os.path.abspath(os.path.expanduser(filename))
1557            self.filename = filename
1558            self.data, self.names = make_empty_project(author=author, filename=filename)
1559        elif os.path.exists(os.path.expanduser(gpxfile)):
1560            # TODO set author
1561            self.data, self.names = LoadDictFromProjFile(gpxfile)
1562            self.update_ids()
1563            if filename:
1564                filename = os.path.abspath(os.path.expanduser(filename))
1565                dr = os.path.split(filename)[0]
1566                if not os.path.exists(dr):
1567                    raise Exception("Directory {} for filename/newgpx does not exist".format(dr))
1568                self.filename = filename
1569            else: 
1570                self.filename = os.path.abspath(os.path.expanduser(gpxfile))
1571        else:
1572            raise ValueError("Not sure what to do with gpxfile {}. Does not exist?".format(gpxfile))
1573
1574    @classmethod
1575    def from_dict_and_names(cls, gpxdict, names, filename=None):
1576        """Creates a :class:`G2Project` directly from
1577        a dictionary and a list of names. If in doubt, do not use this.
1578
1579        :returns: a :class:`G2Project`
1580        """
1581        out = cls()
1582        if filename:
1583            filename = os.path.abspath(os.path.expanduser(filename))
1584            out.filename = filename
1585            gpxdict['Controls']['data']['LastSavedAs'] = filename
1586        else:
1587            try:
1588                out.filename = gpxdict['Controls']['data']['LastSavedAs']
1589            except KeyError:
1590                out.filename = None
1591        out.data = gpxdict
1592        out.names = names
1593
1594    def save(self, filename=None):
1595        """Saves the project, either to the current filename, or to a new file.
1596
1597        Updates self.filename if a new filename provided"""
1598        # TODO update LastSavedUsing ?
1599        if filename:
1600            filename = os.path.abspath(os.path.expanduser(filename))
1601            self.data['Controls']['data']['LastSavedAs'] = filename
1602            self.filename = filename
1603        elif not self.filename:
1604            raise AttributeError("No file name to save to")
1605        SaveDictToProjFile(self.data, self.names, self.filename)
1606
1607    def add_powder_histogram(self, datafile, iparams, phases=[], fmthint=None,
1608                                 databank=None, instbank=None):
1609        """Loads a powder data histogram into the project.
1610
1611        Automatically checks for an instrument parameter file, or one can be
1612        provided. Note that in unix fashion, "~" can be used to indicate the
1613        home directory (e.g. ~/G2data/data.fxye).
1614
1615        :param str datafile: The powder data file to read, a filename.
1616        :param str iparams: The instrument parameters file, a filename.
1617        :param list phases: A list of phases to link to the new histogram,
1618           phases can be references by object, name, rId or number.
1619           Alternately, use 'all' to link to all phases in the project.
1620        :param str fmthint: If specified, only importers where the format name
1621          (reader.formatName, as shown in Import menu) contains the
1622          supplied string will be tried as importers. If not specified, all
1623          importers consistent with the file extension will be tried
1624          (equivalent to "guess format" in menu).
1625        :param int databank: Specifies a dataset number to read, if file contains
1626          more than set of data. This should be 1 to read the first bank in
1627          the file (etc.) regardless of the number on the Bank line, etc.
1628          Default is None which means there should only be one dataset in the
1629          file.
1630        :param int instbank: Specifies an instrument parameter set to read, if
1631          the instrument parameter file contains more than set of parameters.
1632          This will match the INS # in an GSAS type file so it will typically
1633          be 1 to read the first parameter set in the file (etc.)
1634          Default is None which means there should only be one parameter set
1635          in the file.
1636
1637        :returns: A :class:`G2PwdrData` object representing
1638            the histogram
1639        """
1640        LoadG2fil()
1641        datafile = os.path.abspath(os.path.expanduser(datafile))
1642        iparams = os.path.abspath(os.path.expanduser(iparams))
1643        pwdrreaders = import_generic(datafile, Readers['Pwdr'],fmthint=fmthint,bank=databank)
1644        histname, new_names, pwdrdata = load_pwd_from_reader(
1645                                          pwdrreaders[0], iparams,
1646                                          [h.name for h in self.histograms()],bank=instbank)
1647        if histname in self.data:
1648            G2fil.G2Print("Warning - redefining histogram", histname)
1649        elif self.names[-1][0] == 'Phases':
1650            self.names.insert(-1, new_names)
1651        else:
1652            self.names.append(new_names)
1653        self.data[histname] = pwdrdata
1654        self.update_ids()
1655
1656        if phases == 'all':
1657            phases = self.phases()
1658        for phase in phases:
1659            phase = self.phase(phase)
1660            self.link_histogram_phase(histname, phase)
1661
1662        return self.histogram(histname)
1663
1664    def add_simulated_powder_histogram(self, histname, iparams, Tmin, Tmax, Tstep,
1665                                       wavelength=None, scale=None, phases=[]):
1666        """Loads a powder data histogram into the project.
1667
1668        Requires an instrument parameter file.
1669        Note that in unix fashion, "~" can be used to indicate the
1670        home directory (e.g. ~/G2data/data.prm). The instrument parameter file
1671        will determine if the histogram is x-ray, CW neutron, TOF, etc. as well
1672        as the instrument type.
1673
1674        :param str histname: A name for the histogram to be created.
1675        :param str iparams: The instrument parameters file, a filename.
1676        :param float Tmin: Minimum 2theta or TOF (ms) for dataset to be simulated
1677        :param float Tmax: Maximum 2theta or TOF (ms) for dataset to be simulated
1678        :param float Tstep: Step size in 2theta or TOF (ms) for dataset to be simulated       
1679        :param float wavelength: Wavelength for CW instruments, overriding the value
1680           in the instrument parameters file if specified.
1681        :param float scale: Histogram scale factor which multiplies the pattern. Note that
1682           simulated noise is added to the pattern, so that if the maximum intensity is
1683           small, the noise will mask the computed pattern. The scale
1684           needs to be a large number for CW neutrons.
1685           The default, None, provides a scale of 1 for x-rays and TOF; 10,000 for CW neutrons.
1686        :param list phases: Phases to link to the new histogram. Use proj.phases() to link to
1687           all defined phases.
1688
1689        :returns: A :class:`G2PwdrData` object representing the histogram
1690        """
1691        LoadG2fil()
1692        iparams = os.path.abspath(os.path.expanduser(iparams))
1693        if not os.path.exists(iparams):
1694            raise G2ScriptException("File does not exist:"+iparams)
1695        rd = G2obj.ImportPowderData( # Initialize a base class reader
1696            extensionlist=tuple(),
1697            strictExtension=False,
1698            formatName = 'Simulate dataset',
1699            longFormatName = 'Compute a simulated pattern')
1700        rd.powderentry[0] = '' # no filename
1701        rd.powderentry[2] = 1 # only one bank
1702        rd.comments.append('This is a dummy dataset for powder pattern simulation')
1703        rd.idstring = histname
1704        #Iparm1, Iparm2 = load_iprms(iparams, rd)
1705        if Tmax < Tmin:
1706            Tmin,Tmax = Tmax,Tmin
1707        Tstep = abs(Tstep)
1708        if 'TOF' in rd.idstring:
1709                N = (np.log(Tmax)-np.log(Tmin))/Tstep
1710                x = np.exp((np.arange(0,N))*Tstep+np.log(Tmin*1000.))
1711                N = len(x)
1712        else:           
1713                N = int((Tmax-Tmin)/Tstep)+1
1714                x = np.linspace(Tmin,Tmax,N,True)
1715                N = len(x)
1716        if N < 3:
1717            raise G2ScriptException("Error: Range is too small or step is too large, <3 points")
1718        rd.powderdata = [
1719            np.array(x), # x-axis values
1720            np.zeros_like(x), # powder pattern intensities
1721            np.ones_like(x), # 1/sig(intensity)^2 values (weights)
1722            np.zeros_like(x), # calc. intensities (zero)
1723            np.zeros_like(x), # calc. background (zero)
1724            np.zeros_like(x), # obs-calc profiles
1725            ]
1726        Tmin = rd.powderdata[0][0]
1727        Tmax = rd.powderdata[0][-1]
1728        histname, new_names, pwdrdata = load_pwd_from_reader(rd, iparams,
1729                                                            [h.name for h in self.histograms()])
1730        if histname in self.data:
1731            G2fil.G2Print("Warning - redefining histogram", histname)
1732        elif self.names[-1][0] == 'Phases':
1733            self.names.insert(-1, new_names)
1734        else:
1735            self.names.append(new_names)
1736        if scale is not None:
1737            pwdrdata['Sample Parameters']['Scale'][0] = scale
1738        elif pwdrdata['Instrument Parameters'][0]['Type'][0].startswith('PNC'):
1739            pwdrdata['Sample Parameters']['Scale'][0] = 10000.
1740        self.data[histname] = pwdrdata
1741        self.update_ids()
1742
1743        for phase in phases:
1744            phase = self.phase(phase)
1745            self.link_histogram_phase(histname, phase)
1746
1747        return self.histogram(histname)
1748   
1749    def add_phase(self, phasefile, phasename=None, histograms=[], fmthint=None):
1750        """Loads a phase into the project from a .cif file
1751
1752        :param str phasefile: The CIF file from which to import the phase.
1753        :param str phasename: The name of the new phase, or None for the default
1754        :param list histograms: The names of the histograms to associate with
1755            this phase. Use proj.histograms() to add to all histograms.
1756        :param str fmthint: If specified, only importers where the format name
1757          (reader.formatName, as shown in Import menu) contains the
1758          supplied string will be tried as importers. If not specified, all
1759          importers consistent with the file extension will be tried
1760          (equivalent to "guess format" in menu).
1761
1762        :returns: A :class:`G2Phase` object representing the
1763            new phase.
1764        """
1765        LoadG2fil()
1766        histograms = [self.histogram(h).name for h in histograms]
1767        phasefile = os.path.abspath(os.path.expanduser(phasefile))
1768
1769        # TODO handle multiple phases in a file
1770        phasereaders = import_generic(phasefile, Readers['Phase'], fmthint=fmthint)
1771        phasereader = phasereaders[0]
1772       
1773        phasename = phasename or phasereader.Phase['General']['Name']
1774        phaseNameList = [p.name for p in self.phases()]
1775        phasename = G2obj.MakeUniqueLabel(phasename, phaseNameList)
1776        phasereader.Phase['General']['Name'] = phasename
1777
1778        if 'Phases' not in self.data:
1779            self.data[u'Phases'] = { 'data': None }
1780        assert phasename not in self.data['Phases'], "phase names should be unique"
1781        self.data['Phases'][phasename] = phasereader.Phase
1782
1783        if phasereader.Constraints:
1784            Constraints = self.data['Constraints']
1785            for i in phasereader.Constraints:
1786                if isinstance(i, dict):
1787                    if '_Explain' not in Constraints:
1788                        Constraints['_Explain'] = {}
1789                    Constraints['_Explain'].update(i)
1790                else:
1791                    Constraints['Phase'].append(i)
1792
1793        data = self.data['Phases'][phasename]
1794        generalData = data['General']
1795        SGData = generalData['SGData']
1796        NShkl = len(G2spc.MustrainNames(SGData))
1797        NDij = len(G2spc.HStrainNames(SGData))
1798        Super = generalData.get('Super', 0)
1799        if Super:
1800            SuperVec = np.array(generalData['SuperVec'][0])
1801        else:
1802            SuperVec = []
1803        UseList = data['Histograms']
1804
1805        for hist in histograms:
1806            self.link_histogram_phase(hist, phasename)
1807
1808        for obj in self.names:
1809            if obj[0] == 'Phases':
1810                phasenames = obj
1811                break
1812        else:
1813            phasenames = [u'Phases']
1814            self.names.append(phasenames)
1815        phasenames.append(phasename)
1816
1817        # TODO should it be self.filename, not phasefile?
1818        SetupGeneral(data, os.path.dirname(phasefile))
1819        self.index_ids()
1820
1821        self.update_ids()
1822        return self.phase(phasename)
1823
1824    def link_histogram_phase(self, histogram, phase):
1825        """Associates a given histogram and phase.
1826
1827        .. seealso::
1828
1829            :meth:`~G2Project.histogram`
1830            :meth:`~G2Project.phase`"""
1831        hist = self.histogram(histogram)
1832        phase = self.phase(phase)
1833
1834        generalData = phase['General']
1835
1836        if hist.name.startswith('HKLF '):
1837            raise NotImplementedError("HKLF not yet supported")
1838        elif hist.name.startswith('PWDR '):
1839            hist['Reflection Lists'][generalData['Name']] = {}
1840            UseList = phase['Histograms']
1841            SGData = generalData['SGData']
1842            NShkl = len(G2spc.MustrainNames(SGData))
1843            NDij = len(G2spc.HStrainNames(SGData))
1844            UseList[hist.name] = SetDefaultDData('PWDR', hist.name, NShkl=NShkl, NDij=NDij)
1845            UseList[hist.name]['hId'] = hist.id
1846            for key, val in [('Use', True), ('LeBail', False),
1847                             ('newLeBail', True),
1848                             ('Babinet', {'BabA': [0.0, False],
1849                                          'BabU': [0.0, False]})]:
1850                if key not in UseList[hist.name]:
1851                    UseList[hist.name][key] = val
1852        else:
1853            raise RuntimeError("Unexpected histogram" + hist.name)
1854
1855    def reload(self):
1856        """Reload self from self.filename"""
1857        data, names = LoadDictFromProjFile(self.filename)
1858        self.names = names
1859        # Need to deep copy the new data file data into the current tree,
1860        # so that any existing G2Phase, or G2PwdrData objects will still be
1861        # valid
1862        _deep_copy_into(from_=data, into=self.data)
1863
1864    def refine(self, newfile=None, printFile=None, makeBack=False):
1865        '''Invoke a refinement for the project. The project is written to
1866        the currently selected gpx file and then either a single or sequential refinement
1867        is performed depending on the setting of 'Seq Data' in Controls
1868        (set in :meth:`get_Controls`).
1869        '''
1870        seqSetting = self.data['Controls']['data'].get('Seq Data',[])
1871        if not seqSetting:
1872            self.index_ids()    # index_ids will automatically save the project
1873            # TODO: migrate to RefineCore G2strMain does not properly use printFile
1874            # G2strMain.RefineCore(Controls,Histograms,Phases,restraintDict,rigidbodyDict,parmDict,varyList,
1875            #      calcControls,pawleyLookup,ifPrint,printFile,dlg)
1876            G2strMain.Refine(self.filename, makeBack=makeBack)
1877        else:
1878            self._seqrefine()
1879        self.reload() # get file from GPX
1880
1881    def _seqrefine(self):
1882        '''Perform a sequential refinement.
1883        '''
1884
1885        self.data['Controls']['data']['ShowCell'] = True
1886        # add to tree item to project, if not present
1887        if 'Sequential results' not in self.data:
1888            self.data['Sequential results'] = {'data':{}}
1889            self.names.append(['Sequential results'])
1890        self.index_ids()    # index_ids will automatically save the project
1891        #GSASIIpath.IPyBreak_base()
1892
1893        # check that constraints are OK
1894        errmsg, warnmsg = G2strIO.ReadCheckConstraints(self.filename)
1895        if errmsg:
1896            G2fil.G2Print('Refinement error',errmsg)
1897            raise Exception('Constraint error')
1898        if warnmsg:
1899            G2fil.G2Print(u'Warning: Conflict between refinment flag settings and constraints:\n'+
1900                  warnmsg+u'\nRefinement not possible')
1901            raise Exception('Constraint error')
1902        OK,Msg = G2strMain.SeqRefine(self.filename,None)
1903
1904    def histogram(self, histname):
1905        """Returns the histogram named histname, or None if it does not exist.
1906
1907        :param histname: The name of the histogram (str), or ranId or index.
1908        :returns: A :class:`G2PwdrData` object, or None if
1909            the histogram does not exist
1910
1911        .. seealso::
1912            :meth:`~G2Project.histograms`
1913            :meth:`~G2Project.phase`
1914            :meth:`~G2Project.phases`
1915        """
1916        if isinstance(histname, G2PwdrData):
1917            if histname.proj == self:
1918                return histname
1919            else:
1920                raise Exception('Histogram object (G2PwdrData) is not in current project')
1921        if histname in self.data:
1922            return G2PwdrData(self.data[histname], self, histname)
1923        try:
1924            # see if histname is an id or ranId
1925            histname = int(histname)
1926        except ValueError:
1927            return
1928
1929        for histogram in self.histograms():
1930            if histogram.id == histname or histogram.ranId == histname:
1931                return histogram
1932
1933    def histograms(self, typ=None):
1934        """Return a list of all histograms, as :class:`G2PwdrData` objects
1935
1936        For now this only finds Powder/Single Xtal histograms, since that is all that is
1937        currently implemented in this module.
1938
1939        :param ste typ: The prefix (type) the histogram such as 'PWDR '. If None
1940          (the default) all known histograms types are found.
1941        :returns: a list of objects
1942
1943        .. seealso::
1944            :meth:`~G2Project.histogram`
1945            :meth:`~G2Project.phase`
1946            :meth:`~G2Project.phases`
1947        """
1948        output = []
1949        # loop through each tree entry. If it is more than one level (more than one item in the
1950        # list of names). then it must be a histogram, unless labeled Phases or Restraints
1951        if typ is None:
1952            for obj in self.names:
1953                if obj[0].startswith('PWDR ') or obj[0].startswith('HKLF '): 
1954                    output.append(self.histogram(obj[0]))
1955        else:
1956            for obj in self.names:
1957                if len(obj) > 1 and obj[0].startswith(typ): 
1958                    output.append(self.histogram(obj[0]))
1959        return output
1960
1961    def phase(self, phasename):
1962        """
1963        Gives an object representing the specified phase in this project.
1964
1965        :param str phasename: A reference to the desired phase. Either the phase
1966            name (str), the phase's ranId, the phase's index (both int) or
1967            a phase object (:class:`G2Phase`)
1968        :returns: A :class:`G2Phase` object
1969        :raises: KeyError
1970
1971        .. seealso::
1972            :meth:`~G2Project.histograms`
1973            :meth:`~G2Project.phase`
1974            :meth:`~G2Project.phases`
1975            """
1976        if isinstance(phasename, G2Phase):
1977            if phasename.proj == self:
1978                return phasename
1979        phases = self.data['Phases']
1980        if phasename in phases:
1981            return G2Phase(phases[phasename], phasename, self)
1982
1983        try:
1984            # phasename should be phase index or ranId
1985            phasename = int(phasename)
1986        except ValueError:
1987            return
1988
1989        for phase in self.phases():
1990            if phase.id == phasename or phase.ranId == phasename:
1991                return phase
1992
1993    def phases(self):
1994        """
1995        Returns a list of all the phases in the project.
1996
1997        :returns: A list of :class:`G2Phase` objects
1998
1999        .. seealso::
2000            :meth:`~G2Project.histogram`
2001            :meth:`~G2Project.histograms`
2002            :meth:`~G2Project.phase`
2003            """
2004        for obj in self.names:
2005            if obj[0] == 'Phases':
2006                return [self.phase(p) for p in obj[1:]]
2007        return []
2008
2009    def _images(self):
2010        """Returns a list of all the phases in the project.
2011        """
2012        return [i[0] for i in self.names if i[0].startswith('IMG ')]
2013   
2014    def image(self, imageRef):
2015        """
2016        Gives an object representing the specified image in this project.
2017
2018        :param str imageRef: A reference to the desired image. Either the Image
2019          tree name (str), the image's index (int) or
2020          a image object (:class:`G2Image`)
2021        :returns: A :class:`G2Image` object
2022        :raises: KeyError
2023
2024        .. seealso::
2025            :meth:`~G2Project.images`
2026        """
2027        if isinstance(imageRef, G2Image):
2028            if imageRef.proj == self:
2029                return imageRef
2030            else:
2031                raise Exception("Image {} not in current selected project".format(imageRef.name))
2032        if imageRef in self._images():
2033            return G2Image(self.data[imageRef], imageRef, self)
2034
2035        try:
2036            # imageRef should be an index
2037            num = int(imageRef)
2038            imageRef = self._images()[num] 
2039            return G2Image(self.data[imageRef], imageRef, self)
2040        except ValueError:
2041            raise Exception("imageRef {} not an object, name or image index in current selected project"
2042                                .format(imageRef))
2043        except IndexError:
2044            raise Exception("imageRef {} out of range (max={}) in current selected project"
2045                                .format(imageRef,len(self._images())-1))
2046           
2047    def images(self):
2048        """
2049        Returns a list of all the images in the project.
2050
2051        :returns: A list of :class:`G2Image` objects
2052        """
2053        return [G2Image(self.data[i],i,self) for i in self._images()]
2054   
2055    def _pdfs(self):
2056        """Returns a list of all the PDF entries in the project.
2057        """
2058        return [i[0] for i in self.names if i[0].startswith('PDF ')]
2059   
2060    def pdf(self, pdfRef):
2061        """
2062        Gives an object representing the specified PDF entry in this project.
2063
2064        :param pdfRef: A reference to the desired image. Either the PDF
2065          tree name (str), the pdf's index (int) or
2066          a PDF object (:class:`G2PDF`)
2067        :returns: A :class:`G2PDF` object
2068        :raises: KeyError
2069
2070        .. seealso::
2071            :meth:`~G2Project.pdfs`
2072            :class:`~G2PDF`
2073        """
2074        if isinstance(pdfRef, G2PDF):
2075            if pdfRef.proj == self:
2076                return pdfRef
2077            else:
2078                raise Exception("PDF {} not in current selected project".format(pdfRef.name))
2079        if pdfRef in self._pdfs():
2080            return G2PDF(self.data[pdfRef], pdfRef, self)
2081
2082        try:
2083            # pdfRef should be an index
2084            num = int(pdfRef)
2085            pdfRef = self._pdfs()[num] 
2086            return G2PDF(self.data[pdfRef], pdfRef, self)
2087        except ValueError:
2088            raise Exception("pdfRef {} not an object, name or PDF index in current selected project"
2089                                .format(pdfRef))
2090        except IndexError:
2091            raise Exception("pdfRef {} out of range (max={}) in current selected project"
2092                                .format(pdfRef,len(self._images())-1))
2093    def pdfs(self):
2094        """
2095        Returns a list of all the PDFs in the project.
2096
2097        :returns: A list of :class:`G2PDF` objects
2098        """
2099        return [G2PDF(self.data[i],i,self) for i in self._pdfs()]
2100       
2101    def copy_PDF(self, PDFobj, histogram):
2102        '''Creates a PDF entry that can be used to compute a PDF
2103        as a copy of settings in an existing PDF (:class:`G2PDF`)
2104        object.
2105        This places an entry in the project but :meth:`G2PDF.calculate`
2106        must be used to actually perform the PDF computation.
2107
2108        :param str PDFobj: A :class:`G2PDF` object which may be
2109          in a separate project.
2110        :param histogram: A reference to a histogram,
2111          which can be reference by object, name, or number.
2112        :returns: A :class:`G2PDF` object for the PDF entry
2113        '''
2114        LoadG2fil()
2115        PDFname = 'PDF ' + self.histogram(histogram).name[4:]
2116        PDFdict = {'data':None}
2117        for i in 'PDF Controls', 'PDF Peaks':
2118            PDFdict[i] = copy.deepcopy(PDFobj[i])
2119        self.names.append([PDFname]+['PDF Controls', 'PDF Peaks'])
2120        self.data[PDFname] = PDFdict
2121        for i in 'I(Q)','S(Q)','F(Q)','G(R)':
2122            self.data[PDFname]['PDF Controls'][i] = []
2123        G2fil.G2Print('Adding "{}" to project'.format(PDFname))
2124        return G2PDF(self.data[PDFname], PDFname, self)
2125       
2126    def add_PDF(self, prmfile, histogram):
2127        '''Creates a PDF entry that can be used to compute a PDF.
2128        Note that this command places an entry in the project,
2129        but :meth:`G2PDF.calculate` must be used to actually perform
2130        the computation.
2131
2132        :param str datafile: The powder data file to read, a filename.
2133        :param histogram: A reference to a histogram,
2134          which can be reference by object, name, or number.
2135        :returns: A :class:`G2PDF` object for the PDF entry
2136        '''
2137       
2138        LoadG2fil()
2139        PDFname = 'PDF ' + self.histogram(histogram).name[4:]
2140        peaks = {'Limits':[1.,5.],'Background':[2,[0.,-0.2*np.pi],False],'Peaks':[]}
2141        Controls = {
2142            'Sample':{'Name':self.histogram(histogram).name,'Mult':1.0},
2143            'Sample Bkg.':{'Name':'','Mult':-1.0,'Refine':False},
2144            'Container':{'Name':'','Mult':-1.0,'Refine':False},
2145            'Container Bkg.':{'Name':'','Mult':-1.0},
2146            'ElList':{},
2147            'Geometry':'Cylinder','Diam':1.0,'Pack':0.50,'Form Vol':0.0,'Flat Bkg':0,
2148            'DetType':'Area detector','ObliqCoeff':0.2,'Ruland':0.025,'QScaleLim':[20,25],
2149            'Lorch':False,'BackRatio':0.0,'Rmax':100.,'noRing':False,'IofQmin':1.0,'Rmin':1.0,
2150            'I(Q)':[],'S(Q)':[],'F(Q)':[],'G(R)':[]}
2151
2152        fo = open(prmfile,'r')
2153        S = fo.readline()
2154        while S:
2155            if '#' in S:
2156                S = fo.readline()
2157                continue
2158            key,val = S.split(':',1)
2159            try:
2160                Controls[key] = eval(val)
2161            except:
2162                Controls[key] = val.strip()
2163            S = fo.readline()
2164        fo.close()
2165        Controls['Sample']['Name'] = self.histogram(histogram).name
2166        for i in 'Sample Bkg.','Container','Container Bkg.':
2167            Controls[i]['Name'] = ''
2168        PDFdict = {'data':None,'PDF Controls':Controls, 'PDF Peaks':peaks}
2169        self.names.append([PDFname]+['PDF Controls', 'PDF Peaks'])
2170        self.data[PDFname] = PDFdict
2171        G2fil.G2Print('Adding "{}" to project'.format(PDFname))
2172        return G2PDF(self.data[PDFname], PDFname, self)
2173
2174    def seqref(self):
2175        """
2176        Returns a sequential refinement results object, if present
2177
2178        :returns: A :class:`G2SeqRefRes` object or None if not present
2179        """
2180        if 'Sequential results' not in self.data: return
2181        return G2SeqRefRes(self.data['Sequential results']['data'], self)
2182   
2183    def update_ids(self):
2184        """Makes sure all phases and histograms have proper hId and pId"""
2185        # Translated from GetUsedHistogramsAndPhasesfromTree,
2186        #   GSASIIdataGUI.py:4107
2187        for i, h in enumerate(self.histograms()):
2188            h.id = i
2189        for i, p in enumerate(self.phases()):
2190            p.id = i
2191
2192    def do_refinements(self, refinements, histogram='all', phase='all',
2193                       outputnames=None, makeBack=False):
2194        """Conducts one or a series of refinements according to the
2195           input provided in parameter refinements. This is a wrapper
2196           around :meth:`iter_refinements`
2197
2198        :param list refinements: A list of dictionaries specifiying changes to be made to
2199            parameters before refinements are conducted.
2200            See the :ref:`Refinement_recipe` section for how this is defined.
2201        :param str histogram: Name of histogram for refinements to be applied
2202            to, or 'all'; note that this can be overridden for each refinement
2203            step via a "histograms" entry in the dict.
2204        :param str phase: Name of phase for refinements to be applied to, or
2205            'all'; note that this can be overridden for each refinement
2206            step via a "phases" entry in the dict.
2207        :param list outputnames: Provides a list of project (.gpx) file names
2208            to use for each refinement step (specifying None skips the save step).
2209            See :meth:`save`.
2210            Note that this can be overridden using an "output" entry in the dict.
2211        :param bool makeBack: determines if a backup ).bckX.gpx) file is made
2212            before a refinement is performed. The default is False.
2213           
2214        To perform a single refinement without changing any parameters, use this
2215        call:
2216
2217        .. code-block::  python
2218       
2219            my_project.do_refinements([])
2220        """
2221       
2222        for proj in self.iter_refinements(refinements, histogram, phase,
2223                                          outputnames, makeBack):
2224            pass
2225        return self
2226
2227    def iter_refinements(self, refinements, histogram='all', phase='all',
2228                         outputnames=None, makeBack=False):
2229        """Conducts a series of refinements, iteratively. Stops after every
2230        refinement and yields this project, to allow error checking or
2231        logging of intermediate results. Parameter use is the same as for
2232        :meth:`do_refinements` (which calls this method).
2233
2234        >>> def checked_refinements(proj):
2235        ...     for p in proj.iter_refinements(refs):
2236        ...         # Track intermediate results
2237        ...         log(p.histogram('0').residuals)
2238        ...         log(p.phase('0').get_cell())
2239        ...         # Check if parameter diverged, nonsense answer, or whatever
2240        ...         if is_something_wrong(p):
2241        ...             raise Exception("I need a human!")
2242
2243           
2244        """
2245        if outputnames:
2246            if len(refinements) != len(outputnames):
2247                raise ValueError("Should have same number of outputs to"
2248                                 "refinements")
2249        else:
2250            outputnames = [None for r in refinements]
2251
2252        for output, refinedict in zip(outputnames, refinements):
2253            if 'histograms' in refinedict:
2254                hist = refinedict['histograms']
2255            else:
2256                hist = histogram
2257            if 'phases' in refinedict:
2258                ph = refinedict['phases']
2259            else:
2260                ph = phase
2261            if 'output' in refinedict:
2262                output = refinedict['output']
2263            self.set_refinement(refinedict, hist, ph)
2264            # Handle 'once' args - refinements that are disabled after this
2265            # refinement
2266            if 'once' in refinedict:
2267                temp = {'set': refinedict['once']}
2268                self.set_refinement(temp, hist, ph)
2269
2270            if output:
2271                self.save(output)
2272
2273            if 'skip' not in refinedict:
2274                self.refine(makeBack=makeBack)
2275            yield self
2276
2277            # Handle 'once' args - refinements that are disabled after this
2278            # refinement
2279            if 'once' in refinedict:
2280                temp = {'clear': refinedict['once']}
2281                self.set_refinement(temp, hist, ph)
2282            if 'call' in refinedict:
2283                fxn = refinedict['call']
2284                if callable(fxn):
2285                    fxn(*refinedict.get('callargs',[self]))
2286                elif callable(eval(fxn)):
2287                    eval(fxn)(*refinedict.get('callargs',[self]))
2288                else:
2289                    raise G2ScriptException("Error: call value {} is not callable".format(fxn))
2290
2291    def set_refinement(self, refinement, histogram='all', phase='all'):
2292        """Apply specified refinements to a given histogram(s) or phase(s).
2293
2294        :param dict refinement: The refinements to be conducted
2295        :param histogram: Specifies either 'all' (default), a single histogram or
2296          a list of histograms. Histograms may be specified as histogram objects
2297          (see :class:`G2PwdrData`), the histogram name (str) or the index number (int)
2298          of the histogram in the project, numbered starting from 0.
2299          Omitting the parameter or the string 'all' indicates that parameters in
2300          all histograms should be set.
2301        :param phase: Specifies either 'all' (default), a single phase or
2302          a list of phases. Phases may be specified as phase objects
2303          (see :class:`G2Phase`), the phase name (str) or the index number (int)
2304          of the phase in the project, numbered starting from 0.
2305          Omitting the parameter or the string 'all' indicates that parameters in
2306          all phases should be set.
2307
2308        Note that refinement parameters are categorized as one of three types:
2309
2310        1. Histogram parameters
2311        2. Phase parameters
2312        3. Histogram-and-Phase (HAP) parameters
2313       
2314        .. seealso::
2315            :meth:`G2PwdrData.set_refinements`
2316            :meth:`G2PwdrData.clear_refinements`
2317            :meth:`G2Phase.set_refinements`
2318            :meth:`G2Phase.clear_refinements`
2319            :meth:`G2Phase.set_HAP_refinements`
2320            :meth:`G2Phase.clear_HAP_refinements`"""
2321
2322        if histogram == 'all':
2323            hists = self.histograms()
2324        elif isinstance(histogram, list) or isinstance(histogram, tuple):
2325            hists = []
2326            for h in histogram:
2327                if isinstance(h, str) or isinstance(h, int):
2328                    hists.append(self.histogram(h))
2329                else:
2330                    hists.append(h)
2331        elif isinstance(histogram, str) or isinstance(histogram, int):
2332            hists = [self.histogram(histogram)]
2333        else:
2334            hists = [histogram]
2335
2336        if phase == 'all':
2337            phases = self.phases()
2338        elif isinstance(phase, list) or isinstance(phase, tuple):
2339            phases = []
2340            for ph in phase:
2341                if isinstance(ph, str) or isinstance(ph, int):
2342                    phases.append(self.phase(ph))
2343                else:
2344                    phases.append(ph)
2345        elif isinstance(phase, str) or isinstance(phase, int):
2346            phases = [self.phase(phase)]
2347        else:
2348            phases = [phase]
2349
2350        pwdr_set = {}
2351        phase_set = {}
2352        hap_set = {}
2353        for key, val in refinement.get('set', {}).items():
2354            # Apply refinement options
2355            if G2PwdrData.is_valid_refinement_key(key):
2356                pwdr_set[key] = val
2357            elif G2Phase.is_valid_refinement_key(key):
2358                phase_set[key] = val
2359            elif G2Phase.is_valid_HAP_refinement_key(key):
2360                hap_set[key] = val
2361            else:
2362                raise ValueError("Unknown refinement key", key)
2363
2364        for hist in hists:
2365            hist.set_refinements(pwdr_set)
2366        for phase in phases:
2367            phase.set_refinements(phase_set)
2368        for phase in phases:
2369            phase.set_HAP_refinements(hap_set, hists)
2370
2371        pwdr_clear = {}
2372        phase_clear = {}
2373        hap_clear = {}
2374        for key, val in refinement.get('clear', {}).items():
2375            # Clear refinement options
2376            if G2PwdrData.is_valid_refinement_key(key):
2377                pwdr_clear[key] = val
2378            elif G2Phase.is_valid_refinement_key(key):
2379                phase_clear[key] = val
2380            elif G2Phase.is_valid_HAP_refinement_key(key):
2381                hap_set[key] = val
2382            else:
2383                raise ValueError("Unknown refinement key", key)
2384
2385        for hist in hists:
2386            hist.clear_refinements(pwdr_clear)
2387        for phase in phases:
2388            phase.clear_refinements(phase_clear)
2389        for phase in phases:
2390            phase.clear_HAP_refinements(hap_clear, hists)
2391
2392    def index_ids(self):
2393        self.save()
2394        return G2strIO.GetUsedHistogramsAndPhases(self.filename)
2395
2396    def add_constraint_raw(self, cons_scope, constr):
2397        """Adds a constraint of type consType to the project.
2398        cons_scope should be one of "Hist", "Phase", "HAP", or "Global".
2399
2400        WARNING it does not check the constraint is well-constructed"""
2401        constrs = self.data['Constraints']['data']
2402        if 'Global' not in constrs:
2403            constrs['Global'] = []
2404        constrs[cons_scope].append(constr)
2405
2406    def hold_many(self, vars, type):
2407        """Apply holds for all the variables in vars, for constraint of a given type.
2408
2409        type is passed directly to add_constraint_raw as consType
2410
2411        :param list vars: A list of variables to hold. Either :class:`GSASIIobj.G2VarObj` objects,
2412            string variable specifiers, or arguments for :meth:`make_var_obj`
2413        :param str type: A string constraint type specifier. See
2414            :class:`G2Project.add_constraint_raw`
2415
2416        """
2417        for var in vars:
2418            if isinstance(var, str):
2419                var = self.make_var_obj(var)
2420            elif not isinstance(var, G2obj.G2VarObj):
2421                var = self.make_var_obj(*var)
2422            self.add_constraint_raw(type, [[1.0, var], None, None, 'h'])
2423
2424    def make_var_obj(self, phase=None, hist=None, varname=None, atomId=None,
2425                     reloadIdx=True):
2426        """Wrapper to create a G2VarObj. Takes either a string representation ("p:h:name:a")
2427        or individual names of phase, histogram, varname, and atomId.
2428
2429        Automatically converts string phase, hist, or atom names into the ID required
2430        by G2VarObj."""
2431
2432        if reloadIdx:
2433            self.index_ids()
2434
2435        # If string representation, short circuit
2436        if hist is None and varname is None and atomId is None:
2437            if isinstance(phase, str) and ':' in phase:
2438                return G2obj.G2VarObj(phase)
2439
2440        # Get phase index
2441        phaseObj = None
2442        if isinstance(phase, G2Phase):
2443            phaseObj = phase
2444            phase = G2obj.PhaseRanIdLookup[phase.ranId]
2445        elif phase in self.data['Phases']:
2446            phaseObj = self.phase(phase)
2447            phaseRanId = phaseObj.ranId
2448            phase = G2obj.PhaseRanIdLookup[phaseRanId]
2449
2450        # Get histogram index
2451        if isinstance(hist, G2PwdrData):
2452            hist = G2obj.HistRanIdLookup[hist.ranId]
2453        elif hist in self.data:
2454            histRanId = self.histogram(hist).ranId
2455            hist = G2obj.HistRanIdLookup[histRanId]
2456
2457        # Get atom index (if any)
2458        if isinstance(atomId, G2AtomRecord):
2459            atomId = G2obj.AtomRanIdLookup[phase][atomId.ranId]
2460        elif phaseObj:
2461            atomObj = phaseObj.atom(atomId)
2462            if atomObj:
2463                atomRanId = atomObj.ranId
2464                atomId = G2obj.AtomRanIdLookup[phase][atomRanId]
2465
2466        return G2obj.G2VarObj(phase, hist, varname, atomId)
2467
2468    def add_image(self, imagefile, fmthint=None, defaultImage=None):
2469        """Load an image into a project
2470
2471        :param str imagefile: The image file to read, a filename.
2472        :param str fmthint: If specified, only importers where the format name
2473          (reader.formatName, as shown in Import menu) contains the
2474          supplied string will be tried as importers. If not specified, all
2475          importers consistent with the file extension will be tried
2476          (equivalent to "guess format" in menu).
2477        :param str defaultImage: The name of an image to use as a default for
2478          setting parameters for the image file to read.
2479
2480        :returns: a list of G2Image object for the added image(s) [this routine
2481          has not yet been tested with files with multiple images...]
2482        """
2483        LoadG2fil()
2484        imagefile = os.path.abspath(os.path.expanduser(imagefile))
2485        readers = import_generic(imagefile, Readers['Image'], fmthint=fmthint)
2486        objlist = []
2487        for rd in readers:
2488            if rd.SciPy:        #was default read by scipy; needs 1 time fixes
2489                G2fil.G2Print('Warning: Image {} read by generic SciPy import. Image parameters likely wrong'.format(imagefile))
2490                #see this: G2IO.EditImageParms(self,rd.Data,rd.Comments,rd.Image,imagefile)
2491                rd.SciPy = False
2492            rd.readfilename = imagefile
2493            if rd.repeatcount == 1 and not rd.repeat: # skip image number if only one in set
2494                rd.Data['ImageTag'] = None
2495            else:
2496                rd.Data['ImageTag'] = rd.repeatcount
2497            rd.Data['formatName'] = rd.formatName
2498            if rd.sumfile:
2499                rd.readfilename = rd.sumfile
2500            # Load generic metadata, as configured
2501            G2fil.GetColumnMetadata(rd)
2502            # Code from G2IO.LoadImage2Tree(rd.readfilename,self,rd.Comments,rd.Data,rd.Npix,rd.Image)
2503            Imax = np.amax(rd.Image)
2504            ImgNames = [i[0] for i in self.names if i[0].startswith('IMG ')]
2505            TreeLbl = 'IMG '+os.path.basename(imagefile)
2506            ImageTag = rd.Data.get('ImageTag')
2507            if ImageTag:
2508                TreeLbl += ' #'+'%04d'%(ImageTag)
2509                imageInfo = (imagefile,ImageTag)
2510            else:
2511                imageInfo = imagefile
2512            TreeName = G2obj.MakeUniqueLabel(TreeLbl,ImgNames)
2513            # MT dict to contain image info
2514            ImgDict = {}
2515            ImgDict['data'] = [rd.Npix,imageInfo]
2516            ImgDict['Comments'] = rd.Comments
2517            if defaultImage:
2518                if isinstance(defaultImage, G2Image):
2519                    if defaultImage.proj == self:
2520                        defaultImage = self.data[defaultImage.name]['data']
2521                    else:
2522                        raise Exception("Image {} not in current selected project".format(defaultImage.name))
2523                elif defaultImage in self._images():
2524                    defaultImage = self.data[defaultImage]['data']
2525                else:
2526                    defaultImage = None
2527            Data = rd.Data
2528            if defaultImage:
2529                Data = copy.copy(defaultImage)
2530                Data['showLines'] = True
2531                Data['ring'] = []
2532                Data['rings'] = []
2533                Data['cutoff'] = 10.
2534                Data['pixLimit'] = 20
2535                Data['edgemin'] = 100000000
2536                Data['calibdmin'] = 0.5
2537                Data['calibskip'] = 0
2538                Data['ellipses'] = []
2539                Data['calibrant'] = ''
2540                Data['GonioAngles'] = [0.,0.,0.]
2541                Data['DetDepthRef'] = False
2542            else:
2543                Data['type'] = 'PWDR'
2544                Data['color'] = GSASIIpath.GetConfigValue('Contour_color','Paired')
2545                if 'tilt' not in Data:          #defaults if not preset in e.g. Bruker importer
2546                    Data['tilt'] = 0.0
2547                    Data['rotation'] = 0.0
2548                    Data['pixLimit'] = 20
2549                    Data['calibdmin'] = 0.5
2550                    Data['cutoff'] = 10.
2551                Data['showLines'] = False
2552                Data['calibskip'] = 0
2553                Data['ring'] = []
2554                Data['rings'] = []
2555                Data['edgemin'] = 100000000
2556                Data['ellipses'] = []
2557                Data['GonioAngles'] = [0.,0.,0.]
2558                Data['DetDepth'] = 0.
2559                Data['DetDepthRef'] = False
2560                Data['calibrant'] = ''
2561                Data['IOtth'] = [5.0,50.0]
2562                Data['LRazimuth'] = [0.,180.]
2563                Data['azmthOff'] = 0.0
2564                Data['outChannels'] = 2500
2565                Data['outAzimuths'] = 1
2566                Data['centerAzm'] = False
2567                Data['fullIntegrate'] = GSASIIpath.GetConfigValue('fullIntegrate',True)
2568                Data['setRings'] = False
2569                Data['background image'] = ['',-1.0]                           
2570                Data['dark image'] = ['',-1.0]
2571                Data['Flat Bkg'] = 0.0
2572                Data['Oblique'] = [0.5,False]
2573            Data['setDefault'] = False
2574            Data['range'] = [(0,Imax),[0,Imax]]
2575            ImgDict['Image Controls'] = Data
2576            ImgDict['Masks'] = {'Points':[],'Rings':[],'Arcs':[],'Polygons':[],
2577                                'Frames':[],'Thresholds':[(0,Imax),[0,Imax]]}
2578            ImgDict['Stress/Strain']  = {'Type':'True','d-zero':[],'Sample phi':0.0,
2579                                             'Sample z':0.0,'Sample load':0.0}
2580            self.names.append([TreeName]+['Comments','Image Controls','Masks','Stress/Strain'])
2581            self.data[TreeName] = ImgDict
2582            del rd.Image
2583            objlist.append(G2Image(self.data[TreeName], TreeName, self))
2584        return objlist
2585   
2586    def get_Controls(self, control):
2587        '''Return project controls settings
2588
2589        :param str control: the item to be returned. See below for allowed values.
2590        :returns: The value for the control.
2591
2592        Allowed values for parameter control:
2593
2594            * cycles: the maximum number of cycles (returns int)
2595            * sequential: the histograms used for a sequential refinement as a list
2596              of histogram names or an empty list when in non-sequential mode.
2597            * seqCopy: returns True or False. True indicates that results from
2598              each sequential fit are used as the starting point for the next
2599              histogram.
2600            * Anything else returns the value in the Controls dict, if present. An
2601              exception is raised if the control value is not present.
2602
2603        .. seealso::
2604            :meth:`set_Controls`
2605        '''
2606        if control == 'cycles':
2607            return self.data['Controls']['data']['max cyc']
2608        elif control == 'sequential':
2609            return self.data['Controls']['data']['Seq Data']
2610        elif control in self.data['Controls']['data']:
2611            return self.data['Controls']['data'][control]
2612        elif control == 'seqCopy':
2613            return self.data['Controls']['data']['Copy2Next']
2614        else:
2615            G2fil.G2Print('Defined Controls:',self.data['Controls']['data'].keys())
2616            raise Exception('{} is an invalid control value'.format(control))
2617       
2618    def set_Controls(self, control, value):
2619        '''Set project controls
2620
2621        :param str control: the item to be set. See below for allowed values.
2622        :param value: the value to be set.
2623
2624        Allowed values for parameter control:
2625
2626            * cycles: sets the maximum number of cycles (value must be int)
2627            * sequential: sets the histograms to be used for a sequential refinement.
2628              Use an empty list to turn off sequential fitting.
2629              The values in the list may be the name of the histogram (a str), or
2630              a ranId or index (int values), see :meth:`histogram`.
2631            * seqCopy: when True, the results from each sequential fit are used as
2632              the starting point for the next. After each fit is is set to False.
2633              Ignored for non-sequential fits.
2634
2635        .. seealso::
2636            :meth:`get_Controls`
2637        '''
2638        if control == 'cycles':
2639            self.data['Controls']['data']['max cyc'] = int(value)
2640        elif control == 'seqCopy':
2641            self.data['Controls']['data']['Copy2Next'] = bool(value)
2642        elif control == 'sequential':
2643            histlist = []
2644            for i,j in enumerate(value):
2645                h = self.histogram(j)
2646                if h:
2647                    histlist.append(h.name)
2648                else:
2649                    raise Exception('item #{} ({}) is an invalid histogram value'
2650                                        .format(i,j))
2651            self.data['Controls']['data']['Seq Data'] = histlist
2652        else:
2653            raise Exception('{} is an invalid control value'.format(control))
2654       
2655    def copyHistParms(self,sourcehist,targethistlist='all',modelist='all'):
2656        '''Copy histogram information from one histogram to others
2657
2658        :param sourcehist: is a histogram object (:class:`G2PwdrData`) or
2659            a histogram name or the index number of the histogram
2660        :param list targethistlist: a list of histograms where each item in the
2661            list can be a histogram object (:class:`G2PwdrData`),
2662            a histogram name or the index number of the histogram.
2663            if the string 'all' (default value), then all histograms in
2664            the project are used.
2665        :param list modelist: May be a list of sections to copy, which may
2666           include 'Background', 'Instrument Parameters', 'Limits' and
2667           'Sample Parameters' (items may be shortened to uniqueness and
2668           capitalization is ignored, so ['b','i','L','s'] will work.)
2669           The default value, 'all' causes the listed sections to
2670
2671        '''
2672        sections = ('Background','Instrument Parameters','Limits',
2673                        'Sample Parameters')
2674        hist_in = self.histogram(sourcehist)
2675        if not hist_in:
2676            raise Exception('{} is not a valid histogram'.format(sourcehist))
2677        if targethistlist == "all":
2678            targethistlist = self.histograms()
2679        if 'all' in modelist:
2680            copysections = sections
2681        else:
2682            copysections = set()
2683            for s in sections:
2684                for m in modelist:
2685                    if s.lower().startswith(m.lower()):
2686                        copysections.add(s)
2687        for h in targethistlist:
2688            hist_out = self.histogram(h)
2689            if not hist_out:
2690                raise Exception('{} is not a valid histogram'.format(h))
2691            for key in copysections: 
2692                hist_out[key] = copy.deepcopy(hist_in[key])
2693
2694    def get_VaryList(self):
2695        '''Returns a list of the refined variables in the
2696        last refinement cycle
2697
2698        :returns: a list of variables or None if no refinement has been
2699          performed.
2700        '''
2701        try:
2702            return self['Covariance']['data']['varyList']
2703        except:
2704            return
2705
2706    def get_ParmList(self):
2707        '''Returns a list of all the parameters defined in the
2708        last refinement cycle
2709
2710        :returns: a list of parameters or None if no refinement has been
2711          performed.
2712        '''
2713        try:
2714            return list(self['Covariance']['data']['parmDict'].keys())
2715        except:
2716            return
2717       
2718    def get_Variable(self,var):
2719        '''Returns the value and standard uncertainty (esd) for a variable
2720        parameters, as defined in the last refinement cycle
2721
2722        :param str var: a variable name of form '<p>:<h>:<name>', such as
2723          ':0:Scale'
2724        :returns: (value,esd) if the parameter is refined or
2725          (value, None) if the variable is in a constraint or is not
2726          refined or None if the parameter is not found.
2727        '''
2728        if var not in self['Covariance']['data']['parmDict']:
2729            return None
2730        val = self['Covariance']['data']['parmDict'][var]
2731        try:
2732            pos = self['Covariance']['data']['varyList'].index(var)
2733            esd = np.sqrt(self['Covariance']['data']['covMatrix'][pos,pos])
2734            return (val,esd)
2735        except ValueError:
2736            return (val,None)
2737
2738    def get_Covariance(self,varList):
2739        '''Returns the values and covariance matrix for a series of variable
2740        parameters. as defined in the last refinement cycle
2741
2742        :param tuple varList: a list of variable names of form '<p>:<h>:<name>'
2743        :returns: (valueList,CovMatrix) where valueList contains the (n) values
2744          in the same order as varList (also length n) and CovMatrix is a
2745          (n x n) matrix. If any variable name is not found in the varyList
2746          then None is returned.
2747
2748        Use this code, where sig provides standard uncertainties for
2749        parameters and where covArray provides the correlation between
2750        off-diagonal terms::
2751
2752            sig = np.sqrt(np.diag(covMatrix))
2753            xvar = np.outer(sig,np.ones_like(sig))
2754            covArray = np.divide(np.divide(covMatrix,xvar),xvar.T)
2755
2756        '''
2757        missing = [i for i in varList if i not in self['Covariance']['data']['varyList']]
2758        if missing:
2759            G2fil.G2Print('Warning: Variable(s) {} were not found in the varyList'.format(missing))
2760            return None
2761        vals = [self['Covariance']['data']['parmDict'][i] for i in varList]
2762        import GSASIImath as G2mth
2763        cov = G2mth.getVCov(varList,
2764                            self['Covariance']['data']['varyList'],
2765                            self['Covariance']['data']['covMatrix'])
2766        return (vals,cov)
2767
2768class G2AtomRecord(G2ObjectWrapper):
2769    """Wrapper for an atom record. Has convenient accessors via @property.
2770
2771    Available accessors: label, type, refinement_flags, coordinates,
2772        occupancy, ranId, id, adp_flag, uiso
2773
2774    Example:
2775
2776    >>> atom = some_phase.atom("O3")
2777    >>> # We can access the underlying data format
2778    >>> atom.data
2779    ['O3', 'O-2', '', ... ]
2780    >>> # We can also use wrapper accessors
2781    >>> atom.coordinates
2782    (0.33, 0.15, 0.5)
2783    >>> atom.refinement_flags
2784    u'FX'
2785    >>> atom.ranId
2786    4615973324315876477
2787    >>> atom.occupancy
2788    1.0
2789    """
2790    def __init__(self, data, indices, proj):
2791        self.data = data
2792        self.cx, self.ct, self.cs, self.cia = indices
2793        self.proj = proj
2794
2795    @property
2796    def label(self):
2797        '''Get the associated atom's label
2798        '''
2799        return self.data[self.ct-1]
2800
2801    @property
2802    def type(self):
2803        '''Get the associated atom's type
2804        '''
2805        return self.data[self.ct]
2806
2807    @property
2808    def refinement_flags(self):
2809        '''Get or set refinement flags for the associated atom
2810        '''
2811        return self.data[self.ct+1]
2812
2813    @refinement_flags.setter
2814    def refinement_flags(self, other):
2815        # Automatically check it is a valid refinement
2816        for c in other:
2817            if c not in ' FXU':
2818                raise ValueError("Invalid atom refinement: ", other)
2819        self.data[self.ct+1] = other
2820
2821    @property
2822    def coordinates(self):
2823        '''Get the associated atom's coordinates
2824        '''
2825        return tuple(self.data[self.cx:self.cx+3])
2826
2827    @property
2828    def occupancy(self):
2829        '''Get or set the associated atom's occupancy fraction
2830        '''
2831        return self.data[self.cx+3]
2832
2833    @occupancy.setter
2834    def occupancy(self, val):
2835        self.data[self.cx+3] = float(val)
2836
2837    @property
2838    def ranId(self):
2839        '''Get the associated atom's Random Id number
2840        '''
2841        return self.data[self.cia+8]
2842
2843    @property
2844    def adp_flag(self):
2845        '''Get the associated atom's iso/aniso setting, 'I' or 'A'
2846        '''
2847        # Either 'I' or 'A'
2848        return self.data[self.cia]
2849
2850    @property
2851    def uiso(self):
2852        '''Get or set the associated atom's Uiso or Uaniso value(s)
2853        '''
2854        if self.adp_flag == 'I':
2855            return self.data[self.cia+1]
2856        else:
2857            return self.data[self.cia+2:self.cia+8]
2858
2859    @uiso.setter
2860    def uiso(self, value):
2861        if self.adp_flag == 'I':
2862            self.data[self.cia+1] = float(value)
2863        else:
2864            assert len(value) == 6
2865            self.data[self.cia+2:self.cia+8] = [float(v) for v in value]
2866
2867class G2PwdrData(G2ObjectWrapper):
2868    """Wraps a Powder Data Histogram."""
2869    def __init__(self, data, proj, name):
2870        self.data = data
2871        self.proj = proj
2872        self.name = name
2873
2874    @staticmethod
2875    def is_valid_refinement_key(key):
2876        valid_keys = ['Limits', 'Sample Parameters', 'Background',
2877                      'Instrument Parameters']
2878        return key in valid_keys
2879
2880    #@property
2881    #def name(self):
2882    #    return self.data['data'][-1]
2883
2884    @property
2885    def ranId(self):
2886        return self.data['data'][0]['ranId']
2887
2888    @property
2889    def residuals(self):
2890        '''Provides a dictionary with with the R-factors for this histogram.
2891        Includes the weighted and unweighted profile terms (R, Rb, wR, wRb, wRmin)
2892        as well as the Bragg R-values for each phase (ph:H:Rf and ph:H:Rf^2).
2893        '''
2894        data = self.data['data'][0]
2895        return {key: data[key] for key in data
2896                if key in ['R', 'Rb', 'wR', 'wRb', 'wRmin']
2897                   or ':' in key}
2898
2899    @property
2900    def InstrumentParameters(self):
2901        '''Provides a dictionary with with the Instrument Parameters
2902        for this histogram.
2903        '''
2904        return self.data['Instrument Parameters'][0]
2905
2906    @property
2907    def SampleParameters(self):
2908        '''Provides a dictionary with with the Sample Parameters
2909        for this histogram.
2910        '''
2911        return self.data['Sample Parameters']
2912
2913    @property
2914    def Background(self):
2915        '''Provides a list with with the Background parameters
2916        for this histogram.
2917
2918        :returns: list containing a list and dict with background values
2919        '''
2920        return self.data['Background']
2921
2922    def add_back_peak(self,pos,int,sig,gam,refflags=[]):
2923        '''Adds a background peak to the Background parameters
2924       
2925        :param float pos: position of peak, a 2theta or TOF value
2926        :param float int: integrated intensity of background peak, usually large
2927        :param float sig: Gaussian width of background peak, usually large
2928        :param float gam: Lorentzian width of background peak, usually unused (small)
2929        :param list refflags: a list of 1 to 4 boolean refinement flags for
2930            pos,int,sig & gam, respectively (use [0,1] to refine int only).
2931            Defaults to [] which means nothing is refined.
2932        '''
2933        if 'peaksList' not in self.Background[1]:
2934            self.Background[1]['peaksList'] = []
2935        flags = 4*[False]
2936        for i,f in enumerate(refflags):
2937            if i>3: break
2938            flags[i] = bool(f)
2939        bpk = []
2940        for i,j in zip((pos,int,sig,gam),flags):
2941            bpk += [float(i),j]
2942        self.Background[1]['peaksList'].append(bpk)
2943        self.Background[1]['nPeaks'] = len(self.Background[1]['peaksList'])
2944
2945    def del_back_peak(self,peaknum):
2946        '''Removes a background peak from the Background parameters
2947       
2948        :param int peaknum: the number of the peak (starting from 0)
2949        '''
2950        npks = self.Background[1].get('nPeaks',0)
2951        if peaknum >= npks:
2952            raise Exception('peak {} not found in histogram {}'.format(peaknum,self.name))
2953        del self.Background[1]['peaksList'][peaknum]
2954        self.Background[1]['nPeaks'] = len(self.Background[1]['peaksList'])
2955       
2956    def ref_back_peak(self,peaknum,refflags=[]):
2957        '''Sets refinement flag for a background peak
2958       
2959        :param int peaknum: the number of the peak (starting from 0)
2960        :param list refflags: a list of 1 to 4 boolean refinement flags for
2961            pos,int,sig & gam, respectively. If a flag is not specified
2962            it defaults to False (use [0,1] to refine int only).
2963            Defaults to [] which means nothing is refined.
2964        '''
2965        npks = self.Background[1].get('nPeaks',0)
2966        if peaknum >= npks:
2967            raise Exception('peak {} not found in histogram {}'.format(peaknum,self.name))
2968        flags = 4*[False]
2969        for i,f in enumerate(refflags):
2970            if i>3: break
2971            flags[i] = bool(f)
2972        for i,f in enumerate(flags):
2973            self.Background[1]['peaksList'][peaknum][2*i+1] = f
2974                   
2975    @property
2976    def id(self):
2977        self.proj.update_ids()
2978        return self.data['data'][0]['hId']
2979
2980    @id.setter
2981    def id(self, val):
2982        self.data['data'][0]['hId'] = val
2983
2984    def fit_fixed_points(self):
2985        """Attempts to apply a background fit to the fixed points currently specified."""
2986        def SetInstParms(Inst):
2987            dataType = Inst['Type'][0]
2988            insVary = []
2989            insNames = []
2990            insVals = []
2991            for parm in Inst:
2992                insNames.append(parm)
2993                insVals.append(Inst[parm][1])
2994                if parm in ['U','V','W','X','Y','Z','SH/L','I(L2)/I(L1)','alpha',
2995                    'beta-0','beta-1','beta-q','sig-0','sig-1','sig-2','sig-q',] and Inst[parm][2]:
2996                        Inst[parm][2] = False
2997            instDict = dict(zip(insNames, insVals))
2998            instDict['X'] = max(instDict['X'], 0.01)
2999            instDict['Y'] = max(instDict['Y'], 0.01)
3000            if 'SH/L' in instDict:
3001                instDict['SH/L'] = max(instDict['SH/L'], 0.002)
3002            return dataType, instDict, insVary
3003
3004        bgrnd = self.data['Background']
3005
3006        # Need our fixed points in order
3007        bgrnd[1]['FixedPoints'].sort(key=lambda pair: pair[0])
3008        X = [x for x, y in bgrnd[1]['FixedPoints']]
3009        Y = [y for x, y in bgrnd[1]['FixedPoints']]
3010
3011        limits = self.data['Limits'][1]
3012        if X[0] > limits[0]:
3013            X = [limits[0]] + X
3014            Y = [Y[0]] + Y
3015        if X[-1] < limits[1]:
3016            X += [limits[1]]
3017            Y += [Y[-1]]
3018
3019        # Some simple lookups
3020        controls = self.proj['Controls']['data']
3021        inst, inst2 = self.data['Instrument Parameters']
3022        pwddata = self.data['data'][1]
3023
3024        # Construct the data for background fitting
3025        xBeg = np.searchsorted(pwddata[0], limits[0])
3026        xFin = np.searchsorted(pwddata[0], limits[1])
3027        xdata = pwddata[0][xBeg:xFin]
3028        ydata = si.interp1d(X,Y)(ma.getdata(xdata))
3029
3030        W = [1]*len(xdata)
3031        Z = [0]*len(xdata)
3032
3033        dataType, insDict, insVary = SetInstParms(inst)
3034        bakType, bakDict, bakVary = G2pwd.SetBackgroundParms(bgrnd)
3035
3036        # Do the fit
3037        data = np.array([xdata, ydata, W, Z, Z, Z])
3038        G2pwd.DoPeakFit('LSQ', [], bgrnd, limits, inst, inst2, data,
3039                        prevVaryList=bakVary, controls=controls)
3040
3041        # Post-fit
3042        parmDict = {}
3043        bakType, bakDict, bakVary = G2pwd.SetBackgroundParms(bgrnd)
3044        parmDict.update(bakDict)
3045        parmDict.update(insDict)
3046        pwddata[3][xBeg:xFin] *= 0
3047        pwddata[5][xBeg:xFin] *= 0
3048        pwddata[4][xBeg:xFin] = G2pwd.getBackground('', parmDict, bakType, dataType, xdata)[0]
3049
3050        # TODO adjust pwddata? GSASIIpwdGUI.py:1041
3051        # TODO update background
3052        self.proj.save()
3053
3054    def getdata(self,datatype):
3055        '''Provides access to the histogram data of the selected data type
3056
3057        :param str datatype: must be one of the following values (case is ignored)
3058       
3059           * 'X': the 2theta or TOF values for the pattern
3060           * 'Yobs': the observed intensity values
3061           * 'Yweight': the weights for each data point (1/sigma**2)
3062           * 'Ycalc': the computed intensity values
3063           * 'Background': the computed background values
3064           * 'Residual': the difference between Yobs and Ycalc (obs-calc)
3065
3066        :returns: an numpy MaskedArray with data values of the requested type
3067       
3068        '''
3069        enums = ['x', 'yobs', 'yweight', 'ycalc', 'background', 'residual']
3070        if datatype.lower() not in enums:
3071            raise G2ScriptException("Invalid datatype = "+datatype+" must be one of "+str(enums))
3072        return self.data['data'][1][enums.index(datatype.lower())]
3073       
3074    def y_calc(self):
3075        return self.data['data'][1][3]
3076
3077    def Export(self,fileroot,extension):
3078        '''Write the histogram into a file. The path is specified by fileroot and
3079        extension.
3080       
3081        :param str fileroot: name of the file, optionally with a path (extension is
3082           ignored)
3083        :param str extension: includes '.', must match an extension in global
3084           exportersByExtension['powder'] or a Exception is raised.
3085        :returns: name of file that was written
3086        '''
3087        if extension not in exportersByExtension.get('powder',[]):
3088            raise G2ScriptException('No Writer for file type = "'+extension+'"')
3089        fil = os.path.abspath(os.path.splitext(fileroot)[0]+extension)
3090        obj = exportersByExtension['powder'][extension]
3091        obj.SetFromArray(hist=self.data,histname=self.name)
3092        obj.Writer(self.name,fil)
3093        return fil
3094   
3095    def plot(self, Yobs=True, Ycalc=True, Background=True, Residual=True):
3096        try:
3097            import matplotlib.pyplot as plt
3098        except ImportError:
3099            G2fil.G2Print('Warning: Unable to import matplotlib, skipping plot')
3100            return
3101        data = self.data['data'][1]
3102        if Yobs:
3103            plt.plot(data[0], data[1], label='Yobs')
3104        if Ycalc:
3105            plt.plot(data[0], data[3], label='Ycalc')
3106        if Background:
3107            plt.plot(data[0], data[4], label='Background')
3108        if Residual:
3109            plt.plot(data[0], data[5], label="Residual")
3110
3111    def get_wR(self):
3112        """returns the overall weighted profile R factor for a histogram
3113       
3114        :returns: a wR value as a percentage or None if not defined
3115        """
3116        return self['data'][0].get('wR')
3117
3118    def set_refinements(self, refs):
3119        """Sets the histogram refinement parameter 'key' to the specification 'value'.
3120
3121        :param dict refs: A dictionary of the parameters to be set. See the
3122                          :ref:`Histogram_parameters_table` table for a description of
3123                          what these dictionaries should be.
3124
3125        :returns: None
3126
3127        """
3128        do_fit_fixed_points = False
3129        for key, value in refs.items():
3130            if key == 'Limits':
3131                old_limits = self.data['Limits'][1]
3132                new_limits = value
3133                if isinstance(new_limits, dict):
3134                    if 'low' in new_limits:
3135                        old_limits[0] = new_limits['low']
3136                    if 'high' in new_limits:
3137                        old_limits[1] = new_limits['high']
3138                else:
3139                    old_limits[0], old_limits[1] = new_limits
3140            elif key == 'Sample Parameters':
3141                sample = self.data['Sample Parameters']
3142                for sparam in value:
3143                    if sparam not in sample:
3144                        raise ValueError("Unknown refinement parameter, "
3145                                         + str(sparam))
3146                    sample[sparam][1] = True
3147            elif key == 'Background':
3148                bkg, peaks = self.data['Background']
3149
3150                # If True or False, just set the refine parameter
3151                if value in (True, False):
3152                    bkg[1] = value
3153                    return
3154
3155                if 'type' in value:
3156                    bkg[0] = value['type']
3157                if 'refine' in value:
3158                    bkg[1] = value['refine']
3159                if 'no. coeffs' in value:
3160                    cur_coeffs = bkg[2]
3161                    n_coeffs = value['no. coeffs']
3162                    if n_coeffs > cur_coeffs:
3163                        for x in range(n_coeffs - cur_coeffs):
3164                            bkg.append(0.0)
3165                    else:
3166                        for _ in range(cur_coeffs - n_coeffs):
3167                            bkg.pop()
3168                    bkg[2] = n_coeffs
3169                if 'coeffs' in value:
3170                    bkg[3:] = value['coeffs']
3171                if 'FixedPoints' in value:
3172                    peaks['FixedPoints'] = [(float(a), float(b))
3173                                            for a, b in value['FixedPoints']]
3174                if value.get('fit fixed points', False):
3175                    do_fit_fixed_points = True
3176                if 'peaks' in value:
3177                    for i,flags in enumerate(value['peaks']):
3178                        self.ref_back_peak(i,flags)
3179
3180            elif key == 'Instrument Parameters':
3181                instrument, secondary = self.data['Instrument Parameters']
3182                for iparam in value:
3183                    try:
3184                        instrument[iparam][2] = True
3185                    except IndexError:
3186                        raise ValueError("Invalid key:", iparam)
3187            else:
3188                raise ValueError("Unknown key:", key)
3189        # Fit fixed points after the fact - ensure they are after fixed points
3190        # are added
3191        if do_fit_fixed_points:
3192            # Background won't be fit if refinement flag not set
3193            orig = self.data['Background'][0][1]
3194            self.data['Background'][0][1] = True
3195            self.fit_fixed_points()
3196            # Restore the previous value
3197            self.data['Background'][0][1] = orig
3198
3199    def clear_refinements(self, refs):
3200        """Clears the refinement parameter 'key' and its associated value.
3201
3202        :param dict refs: A dictionary of parameters to clear.
3203          See the :ref:`Histogram_parameters_table` table for what can be specified.
3204        """
3205        for key, value in refs.items():
3206            if key == 'Limits':
3207                old_limits, cur_limits = self.data['Limits']
3208                cur_limits[0], cur_limits[1] = old_limits
3209            elif key == 'Sample Parameters':
3210                sample = self.data['Sample Parameters']
3211                for sparam in value:
3212                    sample[sparam][1] = False
3213            elif key == 'Background':
3214                bkg, peaks = self.data['Background']
3215
3216                # If True or False, just set the refine parameter
3217                if value in (True, False):
3218                    bkg[1] = False
3219                    return
3220
3221                bkg[1] = False
3222                if 'FixedPoints' in value:
3223                    if 'FixedPoints' in peaks:
3224                        del peaks['FixedPoints']
3225                if 'peaks' in value:
3226                    for i in range(len(self.Background[1]['peaksList'])):
3227                        self.ref_back_peak(i,[])
3228            elif key == 'Instrument Parameters':
3229                instrument, secondary = self.data['Instrument Parameters']
3230                for iparam in value:
3231                    instrument[iparam][2] = False
3232            else:
3233                raise ValueError("Unknown key:", key)
3234
3235    def add_peak(self,area,dspace=None,Q=None,ttheta=None):
3236        '''Adds a single peak to the peak list
3237        :param float area: peak area
3238        :param float dspace: peak position as d-space (A)
3239        :param float Q: peak position as Q (A-1)
3240        :param float ttheta: peak position as 2Theta (deg)
3241
3242        Note: only one of the parameters dspace, Q or ttheta may be specified
3243        '''
3244        import GSASIIlattice as G2lat
3245        import GSASIImath as G2mth
3246        if (not dspace) + (not Q) + (not ttheta) != 2:
3247            G2fil.G2Print('add_peak error: too many or no peak position(s) specified')
3248            return
3249        pos = ttheta
3250        Parms,Parms2 = self.data['Instrument Parameters']
3251        if Q:
3252            pos = G2lat.Dsp2pos(Parms,2.*np.pi/Q)
3253        elif dspace:
3254            pos = G2lat.Dsp2pos(Parms,dspace)
3255        peaks = self.data['Peak List']
3256        peaks['sigDict'] = {}        #no longer valid
3257        peaks['peaks'].append(G2mth.setPeakparms(Parms,Parms2,pos,area))
3258
3259    def set_peakFlags(self,peaklist=None,area=None,pos=None,sig=None,gam=None):
3260        '''Set refinement flags for peaks
3261       
3262        :param list peaklist: a list of peaks to change flags. If None (default), changes
3263          are made to all peaks.
3264        :param bool area: Sets or clears the refinement flag for the peak area value.
3265          If None (the default), no change is made.
3266        :param bool pos: Sets or clears the refinement flag for the peak position value.
3267          If None (the default), no change is made.
3268        :param bool sig: Sets or clears the refinement flag for the peak sig (Gaussian width) value.
3269          If None (the default), no change is made.
3270        :param bool gam: Sets or clears the refinement flag for the peak sig (Lorentzian width) value.
3271          If None (the default), no change is made.
3272         
3273        Note that when peaks are first created the area flag is on and the other flags are
3274        initially off.
3275
3276        Example::
3277       
3278           set_peakFlags(sig=False,gam=True)
3279
3280        causes the sig refinement flag to be cleared and the gam flag to be set, in both cases for
3281        all peaks. The position and area flags are not changed from their previous values.
3282        '''
3283        peaks = self.data['Peak List']
3284        if peaklist is None:
3285            peaklist = range(len(peaks['peaks']))
3286        for i in peaklist:
3287            for var,j in [(area,3),(pos,1),(sig,5),(gam,7)]:
3288                if var is not None:
3289                    peaks['peaks'][i][j] = var
3290           
3291    def refine_peaks(self):
3292        '''Causes a refinement of peak position, background and instrument parameters
3293        '''
3294        import GSASIIpwd as G2pwd
3295        controls = self.proj.data.get('Controls',{})
3296        controls = controls.get('data',
3297                            {'deriv type':'analytic','min dM/M':0.001,}     #fill in defaults if needed
3298                            )
3299        peaks = self.data['Peak List']
3300        Parms,Parms2 = self.data['Instrument Parameters']
3301        background = self.data['Background']
3302        limits = self.data['Limits'][1]
3303        bxye = np.zeros(len(self.data['data'][1][1]))
3304        peaks['sigDict'] = G2pwd.DoPeakFit('LSQ',peaks['peaks'],background,limits,
3305                                           Parms,Parms2,self.data['data'][1],bxye,[],
3306                                           False,controls,None)[0]
3307
3308    def SaveProfile(self,filename):
3309        '''Writes a GSAS-II (new style) .instprm file
3310        '''
3311        data,Parms2 = self.data['Instrument Parameters']
3312        filename = os.path.splitext(filename)[0]+'.instprm'         # make sure extension is .instprm
3313        File = open(filename,'w')
3314        File.write("#GSAS-II instrument parameter file; do not add/delete items!\n")
3315        for item in data:
3316            File.write(item+':'+str(data[item][1])+'\n')
3317        File.close()
3318        G2fil.G2Print ('Instrument parameters saved to: '+filename)
3319
3320    def LoadProfile(self,filename,bank=0):
3321        '''Reads a GSAS-II (new style) .instprm file and overwrites the current
3322        parameters
3323
3324        :param str filename: instrument parameter file name, extension ignored if not
3325          .instprm
3326        :param int bank: bank number to read, defaults to zero
3327        '''
3328        filename = os.path.splitext(filename)[0]+'.instprm'         # make sure extension is .instprm
3329        File = open(filename,'r')
3330        S = File.readline()
3331        newItems = []
3332        newVals = []
3333        Found = False
3334        while S:
3335            if S[0] == '#':
3336                if Found:
3337                    break
3338                if 'Bank' in S:
3339                    if bank == int(S.split(':')[0].split()[1]):
3340                        S = File.readline()
3341                        continue
3342                    else:
3343                        S = File.readline()
3344                        while S and '#Bank' not in S:
3345                            S = File.readline()
3346                        continue
3347                else:   #a non #Bank file
3348                    S = File.readline()
3349                    continue
3350            Found = True
3351            [item,val] = S[:-1].split(':')
3352            newItems.append(item)
3353            try:
3354                newVals.append(float(val))
3355            except ValueError:
3356                newVals.append(val)                       
3357            S = File.readline()               
3358        File.close()
3359        LoadG2fil()
3360        self.data['Instrument Parameters'][0] = G2fil.makeInstDict(newItems,newVals,len(newVals)*[False,])
3361
3362       
3363class G2Phase(G2ObjectWrapper):
3364    """A wrapper object around a given phase.
3365
3366    Author: Jackson O'Donnell (jacksonhodonnell .at. gmail.com)
3367    """
3368    def __init__(self, data, name, proj):
3369        self.data = data
3370        self.name = name
3371        self.proj = proj
3372
3373    @staticmethod
3374    def is_valid_refinement_key(key):
3375        valid_keys = ["Cell", "Atoms", "LeBail"]
3376        return key in valid_keys
3377
3378    @staticmethod
3379    def is_valid_HAP_refinement_key(key):
3380        valid_keys = ["Babinet", "Extinction", "HStrain", "Mustrain",
3381                      "Pref.Ori.", "Show", "Size", "Use", "Scale"]
3382        return key in valid_keys
3383
3384    def atom(self, atomlabel):
3385        """Returns the atom specified by atomlabel, or None if it does not
3386        exist.
3387
3388        :param str atomlabel: The name of the atom (e.g. "O2")
3389        :returns: A :class:`G2AtomRecord` object
3390            representing the atom.
3391        """
3392        # Consult GSASIIobj.py for the meaning of this
3393        cx, ct, cs, cia = self.data['General']['AtomPtrs']
3394        ptrs = [cx, ct, cs, cia]
3395        atoms = self.data['Atoms']
3396        for atom in atoms:
3397            if atom[ct-1] == atomlabel:
3398                return G2AtomRecord(atom, ptrs, self.proj)
3399
3400    def atoms(self):
3401        """Returns a list of atoms present in the current phase.
3402
3403        :returns: A list of :class:`G2AtomRecord` objects.
3404
3405        .. seealso::
3406            :meth:`~G2Phase.atom`
3407            :class:`G2AtomRecord`
3408        """
3409        ptrs = self.data['General']['AtomPtrs']
3410        output = []
3411        atoms = self.data['Atoms']
3412        for atom in atoms:
3413            output.append(G2AtomRecord(atom, ptrs, self.proj))
3414        return output
3415
3416    def histograms(self):
3417        '''Returns a list of histogram names associated with the current phase ordered
3418        as they appear in the tree (see :meth:`G2Project.histograms`).
3419        '''
3420        return [i.name for i in self.proj.histograms() if i.name in self.data.get('Histograms', {})] 
3421
3422    @property
3423    def ranId(self):
3424        return self.data['ranId']
3425
3426    @property
3427    def id(self):
3428        return self.data['pId']
3429
3430    @id.setter
3431    def id(self, val):
3432        self.data['pId'] = val
3433
3434    def get_cell(self):
3435        """Returns a dictionary of the cell parameters, with keys:
3436            'length_a', 'length_b', 'length_c', 'angle_alpha', 'angle_beta', 'angle_gamma', 'volume'
3437
3438        :returns: a dict
3439
3440        .. seealso::
3441           :meth:`~G2Phase.get_cell_and_esd`
3442
3443        """
3444        cell = self.data['General']['Cell']
3445        return {'length_a': cell[1], 'length_b': cell[2], 'length_c': cell[3],
3446                'angle_alpha': cell[4], 'angle_beta': cell[5], 'angle_gamma': cell[6],
3447                'volume': cell[7]}
3448
3449    def get_cell_and_esd(self):
3450        """
3451        Returns a pair of dictionaries, the first representing the unit cell, the second
3452        representing the estimated standard deviations of the unit cell.
3453
3454        :returns: a tuple of two dictionaries
3455
3456        .. seealso::
3457           :meth:`~G2Phase.get_cell`
3458
3459        """
3460        # translated from GSASIIstrIO.ExportBaseclass.GetCell
3461        import GSASIIlattice as G2lat
3462        import GSASIImapvars as G2mv
3463        try:
3464            pfx = str(self.id) + '::'
3465            sgdata = self['General']['SGData']
3466            covDict = self.proj['Covariance']['data']
3467
3468            parmDict = dict(zip(covDict.get('varyList',[]),
3469                                covDict.get('variables',[])))
3470            sigDict = dict(zip(covDict.get('varyList',[]),
3471                               covDict.get('sig',[])))
3472
3473            if covDict.get('covMatrix') is not None:
3474                sigDict.update(G2mv.ComputeDepESD(covDict['covMatrix'],
3475                                                  covDict['varyList'],
3476                                                  parmDict))
3477
3478            A, sigA = G2strIO.cellFill(pfx, sgdata, parmDict, sigDict)
3479            cellSig = G2strIO.getCellEsd(pfx, sgdata, A, self.proj['Covariance']['data'])
3480            cellList = G2lat.A2cell(A) + (G2lat.calc_V(A),)
3481            cellDict, cellSigDict = {}, {}
3482            for i, key in enumerate(['length_a', 'length_b', 'length_c',
3483                                     'angle_alpha', 'angle_beta', 'angle_gamma',
3484                                     'volume']):
3485                cellDict[key] = cellList[i]
3486                cellSigDict[key] = cellSig[i]
3487            return cellDict, cellSigDict
3488        except KeyError:
3489            cell = self.get_cell()
3490            return cell, {key: 0.0 for key in cell}
3491
3492    def export_CIF(self, outputname, quickmode=True):
3493        """Write this phase to a .cif file named outputname
3494
3495        :param str outputname: The name of the .cif file to write to
3496        :param bool quickmode: Currently ignored. Carryover from exports.G2export_CIF"""
3497        # This code is all taken from exports/G2export_CIF.py
3498        # Functions copied have the same names
3499        import GSASIImath as G2mth
3500        import GSASIImapvars as G2mv
3501        from exports import G2export_CIF as cif
3502
3503        CIFdate = dt.datetime.strftime(dt.datetime.now(),"%Y-%m-%dT%H:%M")
3504        CIFname = os.path.splitext(self.proj.filename)[0]
3505        CIFname = os.path.split(CIFname)[1]
3506        CIFname = ''.join([c if ord(c) < 128 else ''
3507                           for c in CIFname.replace(' ', '_')])
3508        try:
3509            author = self.proj['Controls']['data'].get('Author','').strip()
3510        except KeyError:
3511            pass
3512        oneblock = True
3513
3514        covDict = self.proj['Covariance']['data']
3515        parmDict = dict(zip(covDict.get('varyList',[]),
3516                            covDict.get('variables',[])))
3517        sigDict = dict(zip(covDict.get('varyList',[]),
3518                           covDict.get('sig',[])))
3519
3520        if covDict.get('covMatrix') is not None:
3521            sigDict.update(G2mv.ComputeDepESD(covDict['covMatrix'],
3522                                              covDict['varyList'],
3523                                              parmDict))
3524
3525        with open(outputname, 'w') as fp:
3526            fp.write(' \n' + 70*'#' + '\n')
3527            cif.WriteCIFitem(fp, 'data_' + CIFname)
3528            # from exports.G2export_CIF.WritePhaseInfo
3529            cif.WriteCIFitem(fp, '\n# phase info for '+str(self.name) + ' follows')
3530            cif.WriteCIFitem(fp, '_pd_phase_name', self.name)
3531            # TODO get esds
3532            cellDict = self.get_cell()
3533            defsigL = 3*[-0.00001] + 3*[-0.001] + [-0.01] # significance to use when no sigma
3534            names = ['length_a','length_b','length_c',
3535                     'angle_alpha','angle_beta ','angle_gamma',
3536                     'volume']
3537            for key, val in cellDict.items():
3538                cif.WriteCIFitem(fp, '_cell_' + key, G2mth.ValEsd(val))
3539
3540            cif.WriteCIFitem(fp, '_symmetry_cell_setting',
3541                         self.data['General']['SGData']['SGSys'])
3542
3543            spacegroup = self.data['General']['SGData']['SpGrp'].strip()
3544            # regularize capitalization and remove trailing H/R
3545            spacegroup = spacegroup[0].upper() + spacegroup[1:].lower().rstrip('rh ')
3546            cif.WriteCIFitem(fp, '_symmetry_space_group_name_H-M', spacegroup)
3547
3548            # generate symmetry operations including centering and center of symmetry
3549            SymOpList, offsetList, symOpList, G2oprList, G2opcodes = G2spc.AllOps(
3550                self.data['General']['SGData'])
3551            cif.WriteCIFitem(fp, 'loop_\n    _space_group_symop_id\n    _space_group_symop_operation_xyz')
3552            for i, op in enumerate(SymOpList,start=1):
3553                cif.WriteCIFitem(fp, '   {:3d}  {:}'.format(i,op.lower()))
3554
3555            # TODO skipped histograms, exports/G2export_CIF.py:880
3556
3557            # report atom params
3558            if self.data['General']['Type'] in ['nuclear','macromolecular']:        #this needs macromolecular variant, etc!
3559                cif.WriteAtomsNuclear(fp, self.data, self.name, parmDict, sigDict, [])
3560                # self._WriteAtomsNuclear(fp, parmDict, sigDict)
3561            else:
3562                raise G2ScriptException("no export for "+str(self.data['General']['Type'])+" coordinates implemented")
3563            # report cell contents
3564            cif.WriteComposition(fp, self.data, self.name, parmDict)
3565            if not quickmode and self.data['General']['Type'] == 'nuclear':      # report distances and angles
3566                # WriteDistances(fp,self.name,SymOpList,offsetList,symOpList,G2oprList)
3567                raise NotImplementedError("only quickmode currently supported")
3568            if 'Map' in self.data['General'] and 'minmax' in self.data['General']['Map']:
3569                cif.WriteCIFitem(fp,'\n# Difference density results')
3570                MinMax = self.data['General']['Map']['minmax']
3571                cif.WriteCIFitem(fp,'_refine_diff_density_max',G2mth.ValEsd(MinMax[0],-0.009))
3572                cif.WriteCIFitem(fp,'_refine_diff_density_min',G2mth.ValEsd(MinMax[1],-0.009))
3573
3574
3575    def set_refinements(self, refs):
3576        """Sets the phase refinement parameter 'key' to the specification 'value'
3577
3578        :param dict refs: A dictionary of the parameters to be set. See the
3579                          :ref:`Phase_parameters_table` table for a description of
3580                          this dictionary.
3581
3582        :returns: None"""
3583        for key, value in refs.items():
3584            if key == "Cell":
3585                self.data['General']['Cell'][0] = value
3586
3587            elif key == "Atoms":
3588                for atomlabel, atomrefinement in value.items():
3589                    if atomlabel == 'all':
3590                        for atom in self.atoms():
3591                            atom.refinement_flags = atomrefinement
3592                    else:
3593                        atom = self.atom(atomlabel)
3594                        if atom is None:
3595                            raise ValueError("No such atom: " + atomlabel)
3596                        atom.refinement_flags = atomrefinement
3597
3598            elif key == "LeBail":
3599                hists = self.data['Histograms']
3600                for hname, hoptions in hists.items():
3601                    if 'LeBail' not in hoptions:
3602                        hoptions['newLeBail'] = bool(True)
3603                    hoptions['LeBail'] = bool(value)
3604            else:
3605                raise ValueError("Unknown key:", key)
3606
3607    def clear_refinements(self, refs):
3608        """Clears a given set of parameters.
3609
3610        :param dict refs: The parameters to clear.
3611          See the :ref:`Phase_parameters_table` table for what can be specified.
3612        """
3613        for key, value in refs.items():
3614            if key == "Cell":
3615                self.data['General']['Cell'][0] = False
3616            elif key == "Atoms":
3617                cx, ct, cs, cia = self.data['General']['AtomPtrs']
3618
3619                for atomlabel in value:
3620                    atom = self.atom(atomlabel)
3621                    # Set refinement to none
3622                    atom.refinement_flags = ' '
3623            elif key == "LeBail":
3624                hists = self.data['Histograms']
3625                for hname, hoptions in hists.items():
3626                    if 'LeBail' not in hoptions:
3627                        hoptions['newLeBail'] = True
3628                    hoptions['LeBail'] = False
3629            else:
3630                raise ValueError("Unknown key:", key)
3631
3632    def set_HAP_refinements(self, refs, histograms='all'):
3633        """Sets the given HAP refinement parameters between the current phase and
3634        the specified histograms.
3635
3636        :param dict refs: A dictionary of the parameters to be set. See
3637                          the :ref:`HAP_parameters_table` table for a description of this
3638                          dictionary.
3639        :param histograms: Either 'all' (default) or a list of the histograms by index, name
3640            or object. The index number is relative to all histograms in the tree, not to
3641            those in the phase.
3642            Histograms not associated with the current phase will be ignored.
3643            whose HAP parameters will be set with this phase. Histogram and phase
3644            must already be associated.
3645        :returns: None
3646        """
3647        if not self.data.get('Histograms',[]):
3648            G2fil.G2Print("Error likely: Phase {} has no linked histograms".format(self.name))
3649            return
3650           
3651        if histograms == 'all':
3652            histograms = self.data['Histograms'].keys()
3653        else:
3654            histograms = [self._decodeHist(h) for h in histograms
3655                          if self._decodeHist(h) in self.data['Histograms']]   
3656        if not histograms:
3657            G2fil.G2Print("Warning: Skipping HAP set for phase {}, no selected histograms".format(self.name))
3658            return
3659        for key, val in refs.items():
3660            if key == 'Babinet':
3661                try:
3662                    sets = list(val)
3663                except ValueError:
3664                    sets = ['BabA', 'BabU']
3665                for param in sets:
3666                    if param not in ['BabA', 'BabU']:
3667                        raise ValueError("Not sure what to do with" + param)
3668                    for h in histograms:
3669                        self.data['Histograms'][h]['Babinet'][param][1] = True
3670            elif key == 'Extinction':
3671                for h in histograms:
3672                    self.data['Histograms'][h]['Extinction'][1] = bool(val)
3673            elif key == 'HStrain':
3674                if isinstance(val,list) or isinstance(val,tuple):
3675                    for h in histograms:
3676                        if len(self.data['Histograms'][h]['HStrain'][1]) != len(val):
3677                            raise Exception('Need {} HStrain terms for phase {} hist {}'
3678                                .format(len(self.data['Histograms'][h]['HStrain'][1]),self.name,h))
3679                        for i,v in enumerate(val):
3680                            self.data['Histograms'][h]['HStrain'][1][i] = bool(v)
3681                else:
3682                    for h in histograms:
3683                        self.data['Histograms'][h]['HStrain'][1] = [bool(val) for p in self.data['Histograms'][h]['HStrain'][1]]
3684            elif key == 'Mustrain':
3685                for h in histograms:
3686                    mustrain = self.data['Histograms'][h]['Mustrain']
3687                    newType = None
3688                    direction = None
3689                    if isinstance(val, strtypes):
3690                        if val in ['isotropic', 'uniaxial', 'generalized']:
3691                            newType = val
3692                        else:
3693                            raise ValueError("Not a Mustrain type: " + val)
3694                    elif isinstance(val, dict):
3695                        newType = val.get('type', None)
3696                        direction = val.get('direction', None)
3697
3698                    if newType:
3699                        mustrain[0] = newType
3700                        if newType == 'isotropic':
3701                            mustrain[2][0] = True == val.get('refine',False)
3702                            mustrain[5] = [False for p in mustrain[4]]
3703                        elif newType == 'uniaxial':
3704                            if 'refine' in val:
3705                                mustrain[2][0] = False
3706                                types = val['refine']
3707                                if isinstance(types, strtypes):
3708                                    types = [types]
3709                                elif isinstance(types, bool):
3710                                    mustrain[2][1] = types
3711                                    mustrain[2][2] = types
3712                                    types = []
3713                                else:
3714                                    raise ValueError("Not sure what to do with: "
3715                                                     + str(types))
3716                            else:
3717                                types = []
3718
3719                            for unitype in types:
3720                                if unitype == 'equatorial':
3721                                    mustrain[2][0] = True
3722                                elif unitype == 'axial':
3723                                    mustrain[2][1] = True
3724                                else:
3725                                    msg = 'Invalid uniaxial mustrain type'
3726                                    raise ValueError(msg + ': ' + unitype)
3727                        else:  # newtype == 'generalized'
3728                            mustrain[2] = [False for p in mustrain[1]]
3729                            if 'refine' in val:
3730                                mustrain[5] = [True == val['refine']]*len(mustrain[5])
3731
3732                    if direction:
3733                        if len(direction) != 3:
3734                            raise ValueError("Expected hkl, found", direction)
3735                        direction = [int(n) for n in direction]
3736                        mustrain[3] = direction
3737            elif key == 'Size':
3738                newSize = None
3739                if 'value' in val:
3740                    newSize = float(val['value'])
3741                for h in histograms:
3742                    size = self.data['Histograms'][h]['Size']
3743                    newType = None
3744                    direction = None
3745                    if isinstance(val, strtypes):
3746                        if val in ['isotropic', 'uniaxial', 'ellipsoidal']:
3747                            newType = val
3748                        else:
3749                            raise ValueError("Not a valid Size type: " + val)
3750                    elif isinstance(val, dict):
3751                        newType = val.get('type', None)
3752                        direction = val.get('direction', None)
3753
3754                    if newType:
3755                        size[0] = newType
3756                        refine = bool(val.get('refine'))
3757                        if newType == 'isotropic' and refine is not None:
3758                            size[2][0] = bool(refine)
3759                            if newSize: size[1][0] = newSize
3760                        elif newType == 'uniaxial' and refine is not None:
3761                            size[2][1] = bool(refine)
3762                            size[2][2] = bool(refine)
3763                            if newSize: size[1][1] = size[1][2] =newSize
3764                        elif newType == 'ellipsoidal' and refine is not None:
3765                            size[5] = [bool(refine) for p in size[5]]
3766                            if newSize: size[4] = [newSize for p in size[4]]
3767
3768                    if direction:
3769                        if len(direction) != 3:
3770                            raise ValueError("Expected hkl, found", direction)
3771                        direction = [int(n) for n in direction]
3772                        size[3] = direction
3773            elif key == 'Pref.Ori.':
3774                for h in histograms:
3775                    self.data['Histograms'][h]['Pref.Ori.'][2] = bool(val)
3776            elif key == 'Show':
3777                for h in histograms:
3778                    self.data['Histograms'][h]['Show'] = bool(val)
3779            elif key == 'Use':
3780                for h in histograms:
3781                    self.data['Histograms'][h]['Use'] = bool(val)
3782            elif key == 'Scale':
3783                for h in histograms:
3784                    self.data['Histograms'][h]['Scale'][1] = bool(val)
3785            else:
3786                G2fil.G2Print(u'Warning: Unknown HAP key: '+key)
3787
3788    def clear_HAP_refinements(self, refs, histograms='all'):
3789        """Clears the given HAP refinement parameters between this phase and
3790        the given histograms.
3791
3792        :param dict refs: A dictionary of the parameters to be cleared.
3793            See the the :ref:`HAP_parameters_table` table for what can be specified.
3794        :param histograms: Either 'all' (default) or a list of the histograms by index, name
3795            or object.
3796            The index number is relative to all histograms in the tree, not to
3797            those in the phase.
3798            Histograms not associated with the current phase will be ignored.
3799            whose HAP parameters will be set with this phase. Histogram and phase
3800            must already be associated
3801        :returns: None
3802        """
3803        if histograms == 'all':
3804            histograms = self.data['Histograms'].keys()
3805        else:
3806            histograms = [self._decodeHist(h) for h in histograms
3807                          if self._decodeHist(h) in self.data['Histograms']]   
3808
3809        for key, val in refs.items():
3810            for h in histograms:
3811                if key == 'Babinet':
3812                    try:
3813                        sets = list(val)
3814                    except ValueError:
3815                        sets = ['BabA', 'BabU']
3816                    for param in sets:
3817                        if param not in ['BabA', 'BabU']:
3818                            raise ValueError("Not sure what to do with" + param)
3819                        for h in histograms:
3820                            self.data['Histograms'][h]['Babinet'][param][1] = False
3821                elif key == 'Extinction':
3822                    for h in histograms:
3823                        self.data['Histograms'][h]['Extinction'][1] = False
3824                elif key == 'HStrain':
3825                    for h in histograms:
3826                        self.data['Histograms'][h]['HStrain'][1] = [False for p in self.data['Histograms'][h]['HStrain'][1]]
3827                elif key == 'Mustrain':
3828                    for h in histograms:
3829                        mustrain = self.data['Histograms'][h]['Mustrain']
3830                        mustrain[2] = [False for p in mustrain[2]]
3831                        mustrain[5] = [False for p in mustrain[4]]
3832                elif key == 'Pref.Ori.':
3833                    for h in histograms:
3834                        self.data['Histograms'][h]['Pref.Ori.'][2] = False
3835                elif key == 'Show':
3836                    for h in histograms:
3837                        self.data['Histograms'][h]['Show'] = False
3838                elif key == 'Size':
3839                    for h in histograms:
3840                        size = self.data['Histograms'][h]['Size']
3841                        size[2] = [False for p in size[2]]
3842                        size[5] = [False for p in size[5]]
3843                elif key == 'Use':
3844                    for h in histograms:
3845                        self.data['Histograms'][h]['Use'] = False
3846                elif key == 'Scale':
3847                    for h in histograms:
3848                        self.data['Histograms'][h]['Scale'][1] = False
3849                else:
3850                    G2fil.G2Print(u'Warning: Unknown HAP key: '+key)
3851
3852    def _decodeHist(self,hist):
3853        '''Convert a histogram reference to a histogram name string
3854        '''
3855        if isinstance(hist, G2PwdrData):
3856            return hist.name
3857        elif hist in self.data['Histograms']:
3858            return hist
3859        elif type(hist) is int:
3860            return self.proj.histograms()[hist].name
3861        else:
3862            raise G2ScriptException("Invalid histogram reference: "+str(hist))
3863       
3864    def getHAPvalues(self, histname):
3865        """Returns a dict with HAP values for the selected histogram
3866
3867        :param histogram: is a histogram object (:class:`G2PwdrData`) or
3868            a histogram name or the index number of the histogram.
3869            The index number is relative to all histograms in the tree, not to
3870            those in the phase.
3871        :returns: HAP value dict
3872        """
3873        return self.data['Histograms'][self._decodeHist(histname)]
3874
3875    def copyHAPvalues(self, sourcehist, targethistlist='all', skip=[], use=None):
3876        """Copies HAP parameters for one histogram to a list of other histograms.
3877        Use skip or use to select specific entries to be copied or not used.
3878
3879        :param sourcehist: is a histogram object (:class:`G2PwdrData`) or
3880            a histogram name or the index number of the histogram to copy
3881            parameters from.
3882            The index number is relative to all histograms in the tree, not to
3883            those in the phase.
3884        :param list targethistlist: a list of histograms where each item in the
3885            list can be a histogram object (:class:`G2PwdrData`),
3886            a histogram name or the index number of the histogram.
3887            If the string 'all' (default), then all histograms in the phase
3888            are used.
3889        :param list skip: items in the HAP dict that should not be
3890            copied. The default is an empty list, which causes all items
3891            to be copied. To see a list of items in the dict, use
3892            :meth:`getHAPvalues` or use an invalid item, such as '?'.
3893        :param list use: specifies the items in the HAP dict should be
3894            copied. The default is None, which causes all items
3895            to be copied.
3896
3897        examples::
3898
3899            ph0.copyHAPvalues(0,[1,2,3])
3900            ph0.copyHAPvalues(0,use=['HStrain','Size'])
3901
3902        The first example copies all HAP parameters from the first histogram to
3903        the second, third and fourth histograms (as listed in the project tree).
3904        The second example copies only the 'HStrain' (Dij parameters and
3905        refinement flags) and the 'Size' (crystallite size settings, parameters
3906        and refinement flags) from the first histogram to all histograms.
3907        """
3908        sourcehist = self._decodeHist(sourcehist)
3909        if targethistlist == 'all':
3910            targethistlist = self.histograms()
3911       
3912        copydict = copy.deepcopy(self.data['Histograms'][sourcehist])
3913        for item in skip:
3914            if item in list(copydict.keys()):
3915                del copydict[item]
3916            else:
3917                G2fil.G2Print('items in HAP dict are: {}'.format(
3918                    list(self.data['Histograms'][sourcehist])))
3919                raise Exception('HAP skip list entry {} invalid'.format(item))
3920        if use:
3921            for item in list(copydict.keys()):
3922                if item not in use:
3923                    del copydict[item]
3924
3925        G2fil.G2Print('Copying item(s) {} from histogram {}'.format(list(copydict.keys()),sourcehist))
3926        G2fil.G2Print(' to histogram(s) {}'.format([self._decodeHist(h) for h in targethistlist]))
3927        for h in targethistlist:
3928            h = self._decodeHist(h)
3929            if h not in self.data['Histograms']:
3930                G2fil.G2Print('Unexpected Warning: histogram {} not in phase {}'.format(h,self.name))
3931                continue
3932            self.data['Histograms'][h].update(copy.deepcopy(copydict))
3933
3934    def setHAPvalues(self, HAPdict, targethistlist='all', skip=[], use=None):
3935        """Copies HAP parameters for one histogram to a list of other histograms.
3936        Use skip or use to select specific entries to be copied or not used.
3937        Note that ``HStrain`` and sometimes ``Mustrain`` values can be specific to
3938        a Laue class and should be copied with care between phases of different
3939        symmetry. A "sanity check" on the number of Dij terms is made if ``HStrain``
3940        values are copied.
3941
3942        :param dict HAPdict: is a dict returned by :meth:`getHAPvalues` containing
3943            HAP parameters.
3944        :param list targethistlist: a list of histograms where each item in the
3945            list can be a histogram object (:class:`G2PwdrData`),
3946            a histogram name or the index number of the histogram.
3947            The index number is relative to all histograms in the tree, not to
3948            those in the phase.
3949            If the string 'all' (default), then all histograms in the phase
3950            are used.
3951        :param list skip: items in the HAP dict that should not be
3952            copied. The default is an empty list, which causes all items
3953            to be copied. To see a list of items in the dict, use
3954            :meth:`getHAPvalues` or use an invalid item, such as '?'.
3955        :param list use: specifies the items in the HAP dict should be
3956            copied. The default is None, which causes all items
3957            to be copied.
3958
3959        example::
3960
3961            HAPdict = ph0.getHAPvalues(0)
3962            ph1.setHAPvalues(HAPdict,use=['HStrain','Size'])
3963
3964        This copies the Dij (hydrostatic strain) HAP parameters and the
3965        crystallite size broadening terms from the first histogram in
3966        phase ``ph0`` to all histograms in phase ``ph1``.
3967        """
3968        if targethistlist == 'all':
3969            targethistlist = self.histograms()
3970        copydict = copy.deepcopy(HAPdict)
3971        for item in skip:
3972            if item in list(copydict.keys()):
3973                del copydict[item]
3974            else:
3975                G2fil.G2Print('items in HAP dict are: {}'.format(
3976                    list(self.data['Histograms'][sourcehist])))
3977                raise Exception('HAP skip list entry {} invalid'.format(item))
3978        if use:
3979            for item in list(copydict.keys()):
3980                if item not in use:
3981                    del copydict[item]
3982
3983        first = True
3984        for h in targethistlist:
3985            h = self._decodeHist(h)
3986            if h not in self.data['Histograms']:
3987                G2fil.G2Print('Warning: histogram {} not in phase {}'.format(h,self.name))
3988                continue
3989            if first:
3990                first = False
3991                if 'HStrain' in self.data['Histograms'][h] and 'HStrain' in copydict:
3992                    if len(copydict['HStrain'][0]) != len(self.data['Histograms'][h]['HStrain'][0]):
3993                        G2fil.G2Print('Error: HStrain has differing numbers of terms. Input: {}, phase {}: {}'.
3994                                  format(len(copydict['HStrain'][0]),
3995                                        self.name,len(self.data['Histograms'][h]['HStrain'][0])))
3996                        raise Exception('HStrain has differing numbers of terms.')
3997            self.data['Histograms'][h].update(copy.deepcopy(copydict))           
3998        G2fil.G2Print('Copied item(s) {} from dict'.format(list(copydict.keys())))
3999        G2fil.G2Print(' to histogram(s) {}'.format([self._decodeHist(h) for h in targethistlist]))
4000
4001class G2SeqRefRes(G2ObjectWrapper):
4002    '''Wrapper for a Sequential Refinement Results tree entry, containing the
4003    results for a refinement
4004
4005    As an example::
4006
4007        from __future__ import division, print_function
4008        import os,sys
4009        sys.path.insert(0,'/Users/toby/software/G2/GSASII')
4010        PathWrap = lambda fil: os.path.join('/Users/toby/Scratch/SeqTut2019Mar',fil)
4011        import GSASIIscriptable as G2sc
4012        gpx = G2sc.G2Project(PathWrap('scr4.gpx'))
4013        seq = gpx.seqref()
4014        lbl = ('a','b','c','alpha','beta','gamma','Volume')
4015        for j,h in enumerate(seq.histograms()):
4016            cell,cellU,uniq = seq.get_cell_and_esd(1,h)
4017            print(h)
4018            print([cell[i] for i in list(uniq)+[6]])
4019            print([cellU[i] for i in list(uniq)+[6]])
4020            print('')
4021        print('printed',[lbl[i] for i in list(uniq)+[6]])
4022
4023    .. seealso::
4024        :meth:`G2Project.seqref`
4025    '''
4026    def __init__(self, data, proj):
4027        self.data = data
4028        self.proj = proj
4029        self.newCellDict = {} # dict with recp. cell tensor & Dij name
4030        # newAtomDict = {} # dict with atom positions; relative & absolute
4031        for name in self.data['histNames']:
4032            self.newCellDict.update(self.data[name].get('newCellDict',{}))
4033            #newAtomDict.update(self.data[name].get('newAtomDict',{}) # dict with atom positions; relative & absolute
4034        #ESDlookup = {self.newCellDict[item][0]:item for item in self.newCellDict}
4035        #Dlookup = {item:self.newCellDict[item][0] for item in self.newCellDict}
4036        # Possible error: the next might need to be data[histNames[0]]['varyList']
4037        #atomLookup = {newAtomDict[item][0]:item for item in newAtomDict if item in self.data['varyList']}
4038        #Dlookup.update({atomLookup[parm]:parm for parm in atomLookup}
4039        #ESDlookup.update({parm:atomLookup[parm] for parm in atomLookup})
4040
4041#    @property
4042    def histograms(self):
4043        '''returns a list of histograms in the squential fit
4044        '''
4045        return self.data['histNames']
4046
4047    def get_cell_and_esd(self,phase,hist):
4048        '''Returns a vector of cell lengths and esd values
4049
4050        :param phase: A phase, which may be specified as a phase object
4051          (see :class:`G2Phase`), the phase name (str) or the index number (int)
4052          of the phase in the project, numbered starting from 0.
4053        :param hist: Specify a histogram or using the histogram name (str)
4054          or the index number (int) of the histogram in the sequential
4055          refinement (not the project), numbered as in in the project tree
4056          starting from 0.
4057        :returns: cell,cellESD,uniqCellIndx where cell (list)
4058          with the unit cell parameters (a,b,c,alpha,beta,gamma,Volume);
4059          cellESD are the standard uncertainties on the 7 unit cell
4060          parameters; and uniqCellIndx is a tuple with indicies for the
4061          unique (non-symmetry determined) unit parameters (e.g.
4062          [0,2] for a,c in a tetragonal cell)
4063        '''
4064        def striphist(var,insChar=''):
4065            'strip a histogram number from a var name'
4066            sv = var.split(':')
4067            if len(sv) <= 1: return var
4068            if sv[1]:
4069                sv[1] = insChar
4070            return ':'.join(sv)
4071        import GSASIIlattice as G2lat
4072        import GSASIIstrIO as G2stIO
4073       
4074        uniqCellLookup = [
4075        [['m3','m3m'],(0,)],
4076        [['3R','3mR'],(0,3)],
4077        [['3','3m1','31m','6/m','6/mmm','4/m','4/mmm'],(0,2)],
4078        [['mmm'],(0,1,2)],
4079        [['2/m'+'a'],(0,1,2,3)],
4080        [['2/m'+'b'],(0,1,2,4)],
4081        [['2/m'+'c'],(0,1,2,5)],
4082        [['-1'],(0,1,2,3,4,5)],
4083        ]
4084
4085        seqData,histData = self.RefData(hist)
4086        hId = histData['data'][0]['hId']
4087        phasedict = self.proj.phase(phase).data
4088        pId = phasedict['pId']
4089        pfx = str(pId)+'::' # prefix for A values from phase
4090        phfx = '%d:%d:'%(pId,hId) # Dij prefix
4091        # get unit cell & symmetry for phase
4092        RecpCellTerms = G2lat.cell2A(phasedict['General']['Cell'][1:7])
4093        zeroDict = {pfx+'A'+str(i):0.0 for i in range(6)}
4094        SGdata = phasedict['General']['SGData']
4095        # determine the cell items not defined by symmetry
4096        laue = SGdata['SGLaue'][:]
4097        if laue == '2/m':
4098            laue += SGdata['SGUniq']
4099        for symlist,celllist in uniqCellLookup:
4100            if laue in symlist:
4101                uniqCellIndx = celllist
4102                break
4103        else: # should not happen
4104            uniqCellIndx = list(range(6))
4105        initialCell = {}
4106        for i in uniqCellIndx:
4107            initialCell[str(pId)+'::A'+str(i)] =  RecpCellTerms[i]
4108
4109        esdLookUp = {}
4110        dLookup = {}
4111        # note that varyList keys are p:h:Dij while newCellDict keys are p::Dij
4112        for nKey in seqData['newCellDict']:
4113            p = nKey.find('::')+1
4114            vKey = nKey[:p] + str(hId) + nKey[p:]
4115            if vKey in seqData['varyList']:
4116                esdLookUp[self.newCellDict[nKey][0]] = nKey
4117                dLookup[nKey] = self.newCellDict[nKey][0]
4118        covData = {'varyList': [dLookup.get(striphist(v),v) for v in seqData['varyList']],
4119                'covMatrix': seqData['covMatrix']}
4120        A = RecpCellTerms[:] # make copy of starting A values
4121
4122        for i,j in enumerate(('D11','D22','D33','D12','D13','D23')):
4123            var = pfx+'A'+str(i)
4124            Dvar = phfx+j
4125            # apply Dij value if non-zero
4126            if Dvar in seqData['parmDict']:
4127                A[i] += seqData['parmDict'][Dvar]
4128            # override with fit result if is Dij varied
4129            try:
4130                A[i] = seqData['newCellDict'][esdLookUp[var]][1] # get refined value
4131            except KeyError:
4132                pass
4133        Albls = [pfx+'A'+str(i) for i in range(6)]
4134        cellDict = dict(zip(Albls,A))
4135
4136       
4137        A,zeros = G2stIO.cellFill(pfx,SGdata,cellDict,zeroDict)
4138        # convert to direct cell
4139        c = G2lat.A2cell(A)
4140        vol = G2lat.calc_V(A)
4141        cE = G2stIO.getCellEsd(pfx,SGdata,A,covData)
4142        return list(c)+[vol],cE,uniqCellIndx
4143   
4144    def get_VaryList(self,hist):
4145        '''Returns a list of the refined variables in the
4146        last refinement cycle for the selected histogram
4147
4148        :param hist: Specify a histogram or using the histogram name (str)
4149          or the index number (int) of the histogram in the sequential
4150          refinement (not the project), numbered starting from 0.
4151        :returns: a list of variables or None if no refinement has been
4152          performed.
4153        '''
4154        try:
4155            seqData,histData = self.RefData(hist)
4156            return seqData['varyList']
4157        except:
4158            return
4159
4160    def get_ParmList(self,hist):
4161        '''Returns a list of all the parameters defined in the
4162        last refinement cycle for the selected histogram
4163
4164        :param hist: Specify a histogram or using the histogram name (str)
4165          or the index number (int) of the histogram in the sequential
4166          refinement (not the project), numbered as in the project tree
4167          starting from 0.
4168        :returns: a list of parameters or None if no refinement has been
4169          performed.
4170        '''
4171        try:
4172            seqData,histData = self.RefData(hist)
4173            return list(seqData['parmDict'].keys())
4174        except:
4175            return
4176
4177    def get_Variable(self,hist,var):
4178        '''Returns the value and standard uncertainty (esd) for a variable
4179        parameters, as defined for the selected histogram
4180        in the last sequential refinement cycle
4181
4182        :param hist: Specify a histogram or using the histogram name (str)
4183          or the index number (int) of the histogram in the sequential
4184          refinement (not the project), numbered as in the project tree
4185          starting from 0.
4186        :param str var: a variable name of form '<p>:<h>:<name>', such as
4187          ':0:Scale'
4188        :returns: (value,esd) if the parameter is refined or
4189          (value, None) if the variable is in a constraint or is not
4190          refined or None if the parameter is not found.
4191        '''
4192        try:
4193            seqData,histData = self.RefData(hist)
4194            val = seqData['parmDict'][var]
4195        except:
4196            return
4197        try:
4198            pos = seqData['varyList'].index(var)
4199            esd = np.sqrt(seqData['covMatrix'][pos,pos])
4200            return (val,esd)
4201        except ValueError:
4202            return (val,None)
4203
4204    def get_Covariance(self,hist,varList):
4205        '''Returns the values and covariance matrix for a series of variable
4206        parameters, as defined for the selected histogram
4207        in the last sequential refinement cycle
4208
4209        :param hist: Specify a histogram or using the histogram name (str)
4210          or the index number (int) of the histogram in the sequential
4211          refinement (not the project), numbered as in the project tree
4212          starting from 0.
4213        :param tuple varList: a list of variable names of form '<p>:<h>:<name>'
4214        :returns: (valueList,CovMatrix) where valueList contains the (n) values
4215          in the same order as varList (also length n) and CovMatrix is a
4216          (n x n) matrix. If any variable name is not found in the varyList
4217          then None is returned.
4218
4219        Use this code, where sig provides standard uncertainties for
4220        parameters and where covArray provides the correlation between
4221        off-diagonal terms::
4222
4223            sig = np.sqrt(np.diag(covMatrix))
4224            xvar = np.outer(sig,np.ones_like(sig))
4225            covArray = np.divide(np.divide(covMatrix,xvar),xvar.T)
4226
4227        '''
4228        try:
4229            seqData,histData = self.RefData(hist)
4230        except:
4231            G2fil.G2Print('Warning: Histogram {} not found in the sequential fit'.format(hist))
4232            return
4233        missing = [i for i in varList if i not in seqData['varyList']]
4234        if missing:
4235            G2fil.G2Print('Warning: Variable(s) {} were not found in the varyList'.format(missing))
4236            return None
4237        vals = [seqData['parmDict'][i] for i in varList]
4238        import GSASIImath as G2mth
4239        cov = G2mth.getVCov(varList,seqData['varyList'],seqData['covMatrix'])
4240        return (vals,cov)
4241   
4242    def RefData(self,hist):
4243        '''Provides access to the output from a particular histogram
4244
4245        :param hist: Specify a histogram or using the histogram name (str)
4246          or the index number (int) of the histogram in the sequential
4247          refinement (not the project), numbered as in the project tree
4248          starting from 0.
4249        :returns: a list of dicts where the first element has sequential
4250          refinement results and the second element has the contents of
4251          the histogram tree items.
4252        '''
4253        try:
4254            hist = self.data['histNames'][hist]
4255        except IndexError:
4256            raise Exception('Histogram #{} is out of range from the Sequential Refinement'
4257                                .format(hist))
4258        except TypeError:
4259            pass
4260        if hist not in self.data['histNames']:
4261            raise Exception('Histogram {} is not included in the Sequential Refinement'
4262                                .format(hist))
4263        return self.data[hist],self.proj.histogram(hist).data
4264
4265class G2PDF(G2ObjectWrapper):
4266    """Wrapper for a PDF tree entry, containing the information needed to
4267    compute a PDF and the S(Q), G(r) etc. after the computation is done.
4268    Note that in a GSASIIscriptable script, instances of G2PDF will be created by
4269    calls to :meth:`G2Project.add_PDF` or :meth:`G2Project.pdf`, not via calls
4270    to :meth:`G2PDF.__init__`.
4271
4272    Example use of :class:`G2PDF`::
4273
4274       gpx.add_PDF('250umSiO2.pdfprm',0)
4275       pdf.set_formula(['Si',1],['O',2])
4276       pdf.set_background('Container',1,-0.21)
4277       for i in range(5):
4278           if pdf.optimize(): break
4279       pdf.calculate()
4280       pdf.export(gpx.filename,'S(Q), pdfGUI')
4281       gpx.save('pdfcalc.gpx')
4282
4283    .. seealso::
4284        :meth:`G2Project.pdf`
4285        :meth:`G2Project.pdfs`
4286    """
4287    def __init__(self, data, name, proj):
4288        self.data = data
4289        self.name = name
4290        self.proj = proj
4291    def set_background(self,btype,histogram,mult=-1.,refine=False):
4292        '''Sets a histogram to be used as the 'Sample Background',
4293        the 'Container' or the 'Container Background.'
4294
4295        :param str btype: Type of background to set, must contain
4296          the string 'samp' for Sample Background', 'cont' and 'back'
4297          for the 'Container Background' or only 'cont' for the
4298          'Container'. Note that capitalization and extra characters
4299          are ignored, so the full strings (such as 'Sample
4300          Background' & 'Container Background') can be used.
4301       
4302        :param histogram: A reference to a histogram,
4303          which can be reference by object, name, or number.
4304        :param float mult: a multiplier for the histogram; defaults
4305          to -1.0
4306        :param bool refine: a flag to enable refinement (only
4307          implemented for 'Sample Background'); defaults to False
4308        '''
4309        if 'samp' in btype.lower():
4310            key = 'Sample Bkg.'
4311        elif 'cont' in btype.lower() and 'back' in btype.lower():
4312            key = 'Container Bkg.'
4313        elif 'cont' in btype.lower():
4314            key = 'Container'
4315        else:
4316            raise Exception('btype = {} is invalid'.format(btype))
4317        self.data['PDF Controls'][key]['Name'] = self.proj.histogram(histogram).name
4318        self.data['PDF Controls'][key]['Mult'] = mult
4319        self.data['PDF Controls'][key]['Refine'] = refine
4320
4321    def set_formula(self,*args):
4322        '''Set the chemical formula for the PDF computation.
4323        Use pdf.set_formula(['Si',1],['O',2]) for SiO2.
4324
4325        :param list item1: The element symbol and number of atoms in formula for first element
4326        :param list item2: The element symbol and number of atoms in formula for second element,...
4327
4328        repeat parameters as needed for all elements in the formula.
4329        '''
4330        powderHist = self.proj.histogram(self.data['PDF Controls']['Sample']['Name'])
4331        inst = powderHist.data['Instrument Parameters'][0]
4332        ElList = self.data['PDF Controls']['ElList']
4333        ElList.clear()
4334        sumVol = 0.
4335        for elem,mult in args:
4336            ElList[elem] = G2elem.GetElInfo(elem,inst)
4337            ElList[elem]['FormulaNo'] = mult
4338            Avol = (4.*np.pi/3.)*ElList[elem]['Drad']**3
4339            sumVol += Avol*ElList[elem]['FormulaNo']
4340        self.data['PDF Controls']['Form Vol'] = max(10.0,sumVol)
4341       
4342    def calculate(self,xydata=None,limits=None,inst=None):
4343        '''Compute the PDF using the current parameters. Results are set
4344        in the PDF object arrays (self.data['PDF Controls']['G(R)'] etc.).
4345        Note that if ``xydata``, is specified, the background histograms(s)
4346        will not be accessed from the project file associated with the current
4347        PDF entry. If ``limits`` and ``inst`` are both specified, no histograms
4348        need be in the current project. However, the self.data['PDF Controls']
4349        sections ('Sample', 'Sample Bkg.','Container Bkg.') must be 
4350        non-blank for the corresponding items to be used from``xydata``.
4351
4352        :param dict xydata: an array containing the Sample's I vs Q, and
4353          any or none of the Sample Background, the Container scattering and
4354          the Container Background. If xydata is None (default), the values are
4355          taken from histograms, as named in the PDF's self.data['PDF Controls']
4356          entries with keys 'Sample', 'Sample Bkg.','Container Bkg.' &
4357          'Container'.
4358        :param list limits: upper and lower Q values to be used for PDF
4359          computation. If None (default), the values are
4360          taken from the Sample histogram's .data['Limits'][1] values.
4361        :param dict inst: The Sample histogram's instrument parameters
4362          to be used for PDF computation. If None (default), the values are
4363          taken from the Sample histogram's .data['Instrument Parameters'][0]
4364          values.
4365        '''
4366        data = self.data['PDF Controls']
4367        if xydata is None:
4368            xydata = {}
4369            for key in 'Sample Bkg.','Container Bkg.','Container','Sample':
4370                name = data[key]['Name'].strip()
4371                if name:
4372                    xydata[key] = self.proj.histogram(name).data['data']
4373        if limits is None:
4374            name = data['Sample']['Name'].strip()
4375            limits = self.proj.histogram(name).data['Limits'][1]
4376        if inst is None:
4377            name = data['Sample']['Name'].strip()
4378            inst = self.proj.histogram(name).data['Instrument Parameters'][0]
4379        G2pwd.CalcPDF(data,inst,limits,xydata)
4380        data['I(Q)'] = xydata['IofQ']
4381        data['S(Q)'] = xydata['SofQ']
4382        data['F(Q)'] = xydata['FofQ']
4383        data['G(R)'] = xydata['GofR']
4384
4385    def optimize(self,showFit=True,maxCycles=5,
4386                     xydata=None,limits=None,inst=None):
4387        '''Optimize the low R portion of G(R) to minimize selected
4388        parameters. Note that this updates the parameters in the settings
4389        (self.data['PDF Controls']) but does not update the PDF object
4390        arrays (self.data['PDF Controls']['G(R)'] etc.) with the computed
4391        values, use :meth:`calculate` after a fit to do that.
4392
4393        :param bool showFit: if True (default) the optimized parameters
4394          are shown before and after the fit, as well as the RMS value
4395          in the minimized region.
4396        :param int maxCycles: the maximum number of least-squares cycles;
4397          defaults to 5.
4398        :returns: the result from the optimizer as True or False, depending
4399          on if the refinement converged.
4400        :param dict xydata: an array containing the Sample's I vs Q, and
4401          any or none of the Sample Background, the Container scattering and
4402          the Container Background. If xydata is None (default), the values are
4403          taken from histograms, as named in the PDF's self.data['PDF Controls']
4404          entries with keys 'Sample', 'Sample Bkg.','Container Bkg.' &
4405          'Container'.
4406        :param list limits: upper and lower Q values to be used for PDF
4407          computation. If None (default), the values are
4408          taken from the Sample histogram's .data['Limits'][1] values.
4409        :param dict inst: The Sample histogram's instrument parameters
4410          to be used for PDF computation. If None (default), the values are
4411          taken from the Sample histogram's .data['Instrument Parameters'][0]
4412          values.
4413        '''
4414        data = self.data['PDF Controls']
4415        if xydata is None:
4416            xydata = {}
4417            for key in 'Sample Bkg.','Container Bkg.','Container','Sample':
4418                name = data[key]['Name'].strip()
4419                if name:
4420                    xydata[key] = self.proj.histogram(name).data['data']
4421        if limits is None:
4422            name = data['Sample']['Name'].strip()
4423            limits = self.proj.histogram(name).data['Limits'][1]
4424        if inst is None:
4425            name = data['Sample']['Name'].strip()
4426            inst = self.proj.histogram(name).data['Instrument Parameters'][0]
4427
4428        res = G2pwd.OptimizePDF(data,xydata,limits,inst,showFit,maxCycles)
4429        return res['success']
4430       
4431    def export(self,fileroot,formats):
4432        '''Write out the PDF-related data (G(r), S(Q),...) into files
4433
4434        :param str fileroot: name of file(s) to be written. The extension
4435          will be ignored and set to .iq, .sq, .fq or .gr depending
4436          on the formats selected.
4437        :param str formats: string specifying the file format(s) to be written,
4438          should contain at least one of the following keywords:
4439          I(Q), S(Q), F(Q), G(r) and/or PDFgui (capitalization and
4440          punctuation is ignored). Note that G(r) and PDFgui should not
4441          be specifed together.
4442        '''
4443        PDFsaves = 5*[False]
4444        PDFentry = self.name
4445        name = self.data['PDF Controls']['Sample']['Name'].strip()
4446        limits = self.proj.histogram(name).data['Limits']
4447        inst = self.proj.histogram(name).data['Instrument Parameters'][0]
4448        for i,lbl in enumerate(['I(Q)', 'S(Q)', 'F(Q)', 'G(r)', 'PDFgui']):
4449            PDFsaves[i] = lbl.lower() in formats.lower()
4450        G2fil.PDFWrite(PDFentry,fileroot,PDFsaves,self.data['PDF Controls'],inst,limits)
4451
4452class G2Image(G2ObjectWrapper):
4453    """Wrapper for an IMG tree entry, containing an image and various metadata.
4454    Note that in a GSASIIscriptable script, instances of G2Image will be created by
4455    calls to :meth:`G2Project.add_image` or :meth:`G2Project.images`, not via calls
4456    to :meth:`G2Image.__init__`.
4457
4458    Example use of G2Image:
4459
4460    >>> gpx = G2sc.G2Project(filename='itest.gpx')
4461    >>> imlst = gpx.add_image(idata,fmthint="TIF")
4462    >>> imlst[0].loadControls('stdSettings.imctrl')
4463    >>> imlst[0].setCalibrant('Si    SRM640c')
4464    >>> imlst[0].loadMasks('stdMasks.immask')
4465    >>> imlst[0].Recalibrate()
4466    >>> imlst[0].setControl('outAzimuths',3)
4467    >>> pwdrList = imlst[0].Integrate()
4468
4469
4470    Sample script using an image:   
4471
4472.. code-block::  python
4473
4474    import os,sys
4475    sys.path.insert(0,'/Users/toby/software/G2/GSASII')
4476    import GSASIIscriptable as G2sc
4477    datadir = "/tmp/V4343_IntegrationError/"
4478    PathWrap = lambda fil: os.path.join(datadir,fil)
4479
4480    gpx = G2sc.G2Project(filename=PathWrap('inttest.gpx'))
4481    imlst = gpx.add_image(PathWrap('Si_free_dc800_1-00000.tif'),fmthint="TIF")
4482    imlst[0].loadControls(PathWrap('Si_free_dc800_1-00000.imctrl'))
4483    pwdrList = imlst[0].Integrate()
4484    gpx.save()
4485
4486    """
4487    # parameters in that can be accessed via setControl. This may need future attention
4488    ControlList = {
4489        'int': ['calibskip', 'pixLimit', 'edgemin', 'outChannels',
4490                    'outAzimuths'],
4491        'float': ['cutoff', 'setdist', 'wavelength', 'Flat Bkg',
4492                      'azmthOff', 'tilt', 'calibdmin', 'rotation',
4493                      'distance', 'DetDepth'],
4494        'bool': ['setRings', 'setDefault', 'centerAzm', 'fullIntegrate',
4495                     'DetDepthRef', 'showLines'],
4496        'str': ['SampleShape', 'binType', 'formatName', 'color',
4497                    'type', ],
4498        'list': ['GonioAngles', 'IOtth', 'LRazimuth', 'Oblique', 'PolaVal',
4499                   'SampleAbs', 'center', 'ellipses', 'linescan',
4500                    'pixelSize', 'range', 'ring', 'rings', 'size', ],
4501        'dict': ['varyList'],
4502        }
4503    '''Defines the items known to exist in the Image Controls tree section
4504    and the item's data types. A few are not included here
4505    ('background image', 'dark image', 'Gain map', and 'calibrant') because
4506    these items have special set routines,
4507    where references to entries are checked to make sure their values are
4508    correct.
4509    ''' 
4510       
4511    def __init__(self, data, name, proj):
4512        self.data = data
4513        self.name = name
4514        self.proj = proj
4515
4516    def setControl(self,arg,value):
4517        '''Set an Image Controls parameter in the current image.
4518        If the parameter is not found an exception is raised.
4519
4520        :param str arg: the name of a parameter (dict entry) in the
4521          image. The parameter must be found in :data:`ControlList`
4522          or an exception is raised.
4523        :param value: the value to set the parameter. The value is
4524          cast as the appropriate type from :data:`ControlList`.
4525        '''
4526        for typ in self.ControlList:
4527            if arg in self.ControlList[typ]: break
4528        else:
4529            G2fil.G2Print('Allowed args:\n',[nam for nam,typ in self.findControl('')])
4530            raise Exception('arg {} not defined in G2Image.setControl'
4531                                .format(arg))
4532        try:
4533            if typ == 'int':
4534                self.data['Image Controls'][arg] = int(value)
4535            elif typ == 'float':
4536                self.data['Image Controls'][arg] = float(value)
4537            elif typ == 'bool':
4538                self.data['Image Controls'][arg] = bool(value)
4539            elif typ == 'str':
4540                self.data['Image Controls'][arg] = str(value)
4541            elif typ == 'list':
4542                self.data['Image Controls'][arg] = list(value)
4543            elif typ == 'dict':
4544                self.data['Image Controls'][arg] = dict(value)
4545            else:
4546                raise Exception('Unknown type {} for arg {} in  G2Image.setControl'
4547                                    .format(typ,arg))
4548        except:
4549            raise Exception('Error formatting value {} as type {} for arg {} in  G2Image.setControl'
4550                                    .format(value,typ,arg))
4551
4552    def getControl(self,arg):
4553        '''Return an Image Controls parameter in the current image.
4554        If the parameter is not found an exception is raised.
4555
4556        :param str arg: the name of a parameter (dict entry) in the
4557          image.
4558        :returns: the value as a int, float, list,...
4559        '''
4560        if arg in self.data['Image Controls']:
4561            return self.data['Image Controls'][arg]
4562        G2fil.G2Print(self.findControl(''))
4563        raise Exception('arg {} not defined in G2Image.getControl'.format(arg))
4564
4565    def findControl(self,arg=''):
4566        '''Finds the Image Controls parameter(s) in the current image
4567        that match the string in arg. Default is '' which returns all
4568        parameters.
4569
4570            Example:
4571
4572            >>> findControl('calib')
4573            [['calibskip', 'int'], ['calibdmin', 'float'], ['calibrant', 'str']]
4574
4575        :param str arg: a string containing part of the name of a
4576          parameter (dict entry) in the image's Image Controls.
4577        :returns: a list of matching entries in form
4578          [['item','type'], ['item','type'],...] where each 'item' string
4579          contains the sting in arg.
4580        '''
4581        matchList = []
4582        for typ in self.ControlList:
4583            for item in self.ControlList[typ]:
4584                if arg in item:
4585                    matchList.append([item,typ])
4586        return matchList
4587
4588    def setCalibrant(self,calib):
4589        '''Set a calibrant for the current image
4590
4591        :param str calib: specifies a calibrant name which must be one of
4592          the entries in file ImageCalibrants.py. This is validated and
4593          an error provides a list of valid choices.
4594        '''
4595        import ImageCalibrants as calFile
4596        if calib in calFile.Calibrants.keys():
4597            self.data['Image Controls']['calibrant'] = calib
4598            return
4599        G2fil.G2Print('Calibrant {} is not valid. Valid calibrants'.format(calib))
4600        for i in calFile.Calibrants.keys():
4601            if i: G2fil.G2Print('\t"{}"'.format(i))
4602       
4603    def setControlFile(self,typ,imageRef,mult=None):
4604        '''Set a image to be used as a background/dark/gain map image
4605
4606        :param str typ: specifies image type, which must be one of:
4607           'background image', 'dark image', 'gain map'; N.B. only the first
4608           four characters must be specified and case is ignored.
4609        :param imageRef: A reference to the desired image. Either the Image
4610          tree name (str), the image's index (int) or
4611          a image object (:class:`G2Image`)
4612        :param float mult: a multiplier to be applied to the image (not used
4613          for 'Gain map'; required for 'background image', 'dark image'
4614        '''
4615        if 'back' in typ.lower():
4616            key = 'background image'
4617            mult = float(mult)
4618        elif 'dark' in typ.lower():
4619            key = 'dark image'
4620            mult = float(mult)
4621        elif 'gain' in typ.lower():
4622            #key = 'Gain map'
4623            if mult is not None:
4624                G2fil.G2Print('Warning: Ignoring multiplier for Gain map')
4625            mult = None
4626        else:
4627            raise Exception("Invalid typ {} for setControlFile".format(typ))
4628        imgNam = self.proj.image(imageRef).name
4629        if mult is None:
4630            self.data['Image Controls']['Gain map'] = imgNam
4631        else:
4632            self.data['Image Controls'][key] = [imgNam,mult]
4633
4634    def loadControls(self,filename=None,imgDict=None):
4635        '''load controls from a .imctrl file
4636
4637        :param str filename: specifies a file to be read, which should end
4638          with .imctrl (defaults to None, meaning parameters are input
4639          with imgDict.)
4640        :param dict imgDict: contains a set of image parameters (defaults to
4641          None, meaning parameters are input with filename.)
4642        '''
4643        if filename:
4644            File = open(filename,'r')
4645            Slines = File.readlines()
4646            File.close()
4647            G2fil.LoadControls(Slines,self.data['Image Controls'])
4648            G2fil.G2Print('file {} read into {}'.format(filename,self.name))
4649        elif imgDict:
4650            self.data['Image Controls'].update(imgDict)
4651            G2fil.G2Print('Image controls set')
4652        else:
4653            raise Exception("loadControls called without imgDict or filename specified")
4654
4655    def saveControls(self,filename):
4656        '''write current controls values to a .imctrl file
4657
4658        :param str filename: specifies a file to write, which should end
4659          with .imctrl
4660        '''
4661        G2fil.WriteControls(filename,self.data['Image Controls'])
4662        G2fil.G2Print('file {} written from {}'.format(filename,self.name))
4663
4664    def getControls(self,clean=False):
4665        '''returns current Image Controls as a dict
4666
4667        :param bool clean: causes the calbration information to be deleted
4668        '''
4669        ImageControls = copy.deepcopy(self.data['Image Controls'])
4670        if clean:
4671            ImageControls['showLines'] = True
4672            ImageControls['ring'] = []
4673            ImageControls['rings'] = []
4674            ImageControls['ellipses'] = []
4675            ImageControls['setDefault'] = False
4676            for i in 'range','size','GonioAngles':
4677                if i in ImageControls: del ImageControls[i]
4678        return ImageControls
4679   
4680    def setControls(self,controlsDict):
4681        '''uses dict from :meth:`getControls` to set Image Controls for current image
4682        '''
4683        self.data['Image Controls'].update(copy.deepcopy(controlsDict))
4684   
4685   
4686    def loadMasks(self,filename,ignoreThreshold=False):
4687        '''load masks from a .immask file
4688
4689        :param str filename: specifies a file to be read, which should end
4690          with .immask
4691        :param bool ignoreThreshold: If True, masks are loaded with
4692          threshold masks. Default is False which means any Thresholds
4693          in the file are ignored.
4694        '''
4695        G2fil.readMasks(filename,self.data['Masks'],ignoreThreshold)
4696        G2fil.G2Print('file {} read into {}'.format(filename,self.name))
4697
4698    def initMasks(self):
4699        '''Initialize Masks, including resetting the Thresholds values
4700        '''
4701        self.data['Masks'] = {'Points':[],'Rings':[],'Arcs':[],'Polygons':[],'Frames':[]}
4702        ImageZ = _getCorrImage(Readers['Image'],self.proj,self)
4703        Imin = max(0.,np.min(ImageZ))
4704        Imax = np.max(ImageZ)
4705        self.data['Masks']['Thresholds'] = [(0,Imax),[Imin,Imax]]
4706       
4707    def getMasks(self):
4708        '''load masks from an IMG tree entry
4709        '''
4710        return self.data['Masks']
4711       
4712    def setMasks(self,maskDict,resetThresholds=False):
4713        '''load masks dict (from :meth:`getMasks`) into current IMG record
4714
4715        :param dict maskDict: specifies a dict with image parameters
4716        :param bool resetThresholds: If True, Threshold Masks in the
4717          dict are ignored. The default is False which means Threshold
4718          Masks are retained.
4719        '''
4720        self.data['Masks'] = copy.deepcopy(maskDict)
4721        if resetThresholds:
4722            ImageZ = _getCorrImage(Readers['Image'],self.proj,self)
4723            Imin = max(0.,np.min(ImageZ))
4724            Imax = np.max(ImageZ)
4725            self.data['Masks']['Thresholds'] [(0,Imax),[Imin,Imax]]
4726       
4727    def getVary(self,*args):
4728        '''Return the refinement flag(s) for Image Controls parameter(s)
4729        in the current image.
4730        If the parameter is not found, an exception is raised.
4731
4732        :param str arg: the name of a refinement parameter in the
4733          varyList for the image. The name should be one of
4734          'dep', 'det-X', 'det-Y', 'dist', 'phi', 'tilt', or 'wave'
4735        :param str arg1: the name of a parameter (dict entry) as before,
4736          optional
4737
4738
4739        :returns: a list of bool value(s)
4740        '''
4741        res = []
4742        for arg in args:
4743            if arg in self.data['Image Controls']['varyList']:
4744                res.append(self.data['Image Controls']['varyList'][arg])
4745            else:
4746                raise Exception('arg {} not defined in G2Image.getVary'.format(arg))
4747        return res
4748   
4749    def setVary(self,arg,value):
4750        '''Set a refinement flag for Image Controls parameter in the
4751        current image.
4752        If the parameter is not found an exception is raised.
4753
4754        :param str arg: the name of a refinement parameter in the
4755          varyList for the image. The name should be one of
4756          'dep', 'det-X', 'det-Y', 'dist', 'phi', 'tilt', or 'wave'
4757        :param str arg: the name of a parameter (dict entry) in the
4758          image. The parameter must be found in :data:`ControlList`
4759          or an exception is raised.
4760        :param value: the value to set the parameter. The value is
4761          cast as the appropriate type from :data:`ControlList`.
4762        '''
4763        if arg in self.data['Image Controls']['varyList']:
4764            self.data['Image Controls']['varyList'][arg] = bool(value)
4765        else:
4766            raise Exception('arg {} not defined in G2Image.setVary'.format(arg))
4767
4768    def Recalibrate(self):
4769        '''Invokes a recalibration fit (same as Image Controls/Calibration/Recalibrate
4770        menu command). Note that for this to work properly, the calibration
4771        coefficients (center, wavelength, ditance & tilts) must be fairly close.
4772        This may produce a better result if run more than once.
4773        '''
4774        LoadG2fil()
4775        ImageZ = _getCorrImage(Readers['Image'],self.proj,self)
4776        G2img.ImageRecalibrate(None,ImageZ,self.data['Image Controls'],self.data['Masks'])
4777
4778    def Integrate(self,name=None):
4779        '''Invokes an image integration (same as Image Controls/Integration/Integrate
4780        menu command). All parameters will have previously been set with Image Controls
4781        so no input is needed here. Note that if integration is performed on an
4782        image more than once, histogram entries may be overwritten. Use the name
4783        parameter to prevent this if desired.
4784
4785        :param str name: base name for created histogram(s). If None (default),
4786          the histogram name is taken from the image name.
4787        :returns: a list of created histogram (:class:`G2PwdrData`) objects.
4788        '''
4789        blkSize = 256   #256 seems to be optimal; will break in polymask if >1024
4790        ImageZ = _getCorrImage(Readers['Image'],self.proj,self)
4791        # do integration
4792        ints,azms,Xvals,cancel = G2img.ImageIntegrate(ImageZ,self.data['Image Controls'],self.data['Masks'],blkSize=blkSize)
4793        # code from here on based on G2IO.SaveIntegration, but places results in the current
4794        # project rather than tree
4795        X = Xvals[:-1]
4796        N = len(X)
4797
4798        data = self.data['Image Controls']
4799        Comments = self.data['Comments']
4800        # make name from image, unless overridden
4801        if name:
4802            if not name.startswith(data['type']+' '):
4803                name = data['type']+' '+name
4804        else:
4805            name = self.name.replace('IMG ',data['type']+' ')
4806        if 'PWDR' in name:
4807            if 'target' in data:
4808                names = ['Type','Lam1','Lam2','I(L2)/I(L1)','Zero','Polariz.','U','V','W','X','Y','Z','SH/L','Azimuth'] 
4809                codes = [0 for i in range(14)]
4810            else:
4811                names = ['Type','Lam','Zero','Polariz.','U','V','W','X','Y','Z','SH/L','Azimuth'] 
4812                codes = [0 for i in range(12)]
4813        elif 'SASD' in name:
4814            names = ['Type','Lam','Zero','Azimuth'] 
4815            codes = [0 for i in range(4)]
4816            X = 4.*np.pi*npsind(X/2.)/data['wavelength']    #convert to q
4817        Xminmax = [X[0],X[-1]]
4818        Azms = []
4819        dazm = 0.
4820        if data['fullIntegrate'] and data['outAzimuths'] == 1:
4821            Azms = [45.0,]                              #a poor man's average?
4822        else:
4823            for i,azm in enumerate(azms[:-1]):
4824                if azm > 360. and azms[i+1] > 360.:
4825                    Azms.append(G2img.meanAzm(azm%360.,azms[i+1]%360.))
4826                else:   
4827                    Azms.append(G2img.meanAzm(azm,azms[i+1]))
4828            dazm = np.min(np.abs(np.diff(azms)))/2.
4829        # pull out integration results and make histograms for each
4830        IntgOutList = []
4831        for i,azm in enumerate(azms[:-1]):
4832            Aname = name+" Azm= %.2f"%((azm+dazm)%360.)
4833            # MT dict to contain histogram
4834            HistDict = {}
4835            histItems = [Aname]
4836            Sample = G2obj.SetDefaultSample()       #set as Debye-Scherrer
4837            Sample['Gonio. radius'] = data['distance']
4838            Sample['Omega'] = data['GonioAngles'][0]
4839            Sample['Chi'] = data['GonioAngles'][1]
4840            Sample['Phi'] = data['GonioAngles'][2]
4841            Sample['Azimuth'] = (azm+dazm)%360.    #put here as bin center
4842            polariz = 0.99    #set default polarization for synchrotron radiation!
4843            for item in Comments:
4844                if 'polariz' in item:
4845                    try:
4846                        polariz = float(item.split('=')[1])
4847                    except:
4848                        polariz = 0.99
4849                for key in ('Temperature','Pressure','Time','FreePrm1','FreePrm2','FreePrm3','Omega',
4850                    'Chi','Phi'):
4851                    if key.lower() in item.lower():
4852                        try:
4853                            Sample[key] = float(item.split('=')[1])
4854                        except:
4855                            pass
4856                # if 'label_prm' in item.lower():
4857                #     for num in ('1','2','3'):
4858                #         if 'label_prm'+num in item.lower():
4859                #             Controls['FreePrm'+num] = item.split('=')[1].strip()
4860            if 'PWDR' in Aname:
4861                if 'target' in data:    #from lab x-ray 2D imaging data
4862                    waves = {'CuKa':[1.54051,1.54433],'TiKa':[2.74841,2.75207],'CrKa':[2.28962,2.29351],
4863                                 'FeKa':[1.93597,1.93991],'CoKa':[1.78892,1.79278],'MoKa':[0.70926,0.713543],
4864                                 'AgKa':[0.559363,0.563775]}
4865                    wave1,wave2 = waves[data['target']]
4866                    parms = ['PXC',wave1,wave2,0.5,0.0,polariz,290.,-40.,30.,6.,-14.,0.0,0.0001,Azms[i]]
4867                else:
4868                    parms = ['PXC',data['wavelength'],0.0,polariz,1.0,-0.10,0.4,0.30,1.0,0.0,0.0001,Azms[i]]
4869            elif 'SASD' in Aname:
4870                Sample['Trans'] = data['SampleAbs'][0]
4871                parms = ['LXC',data['wavelength'],0.0,Azms[i]]
4872            Y = ints[i]
4873            Ymin = np.min(Y)
4874            Ymax = np.max(Y)
4875            W = np.where(Y>0.,1./Y,1.e-6)                    #probably not true
4876            section = 'Comments'
4877            histItems += [section]
4878            HistDict[section] = Comments
4879            section = 'Limits'
4880            histItems += [section]
4881            HistDict[section] = copy.deepcopy([tuple(Xminmax),Xminmax])
4882            if 'PWDR' in Aname:
4883                section = 'Background'
4884                histItems += [section]
4885                HistDict[section] = [['chebyschev',1,3,1.0,0.0,0.0],
4886                    {'nDebye':0,'debyeTerms':[],'nPeaks':0,'peaksList':[]}]
4887            inst = [dict(zip(names,zip(parms,parms,codes))),{}]
4888            for item in inst[0]:
4889                inst[0][item] = list(inst[0][item])
4890            section = 'Instrument Parameters'
4891            histItems += [section]
4892            HistDict[section] = inst
4893            if 'PWDR' in Aname:
4894                section = 'Sample Parameters'
4895                histItems += [section]
4896                HistDict[section] = Sample
4897                section = 'Peak List'
4898                histItems += [section]
4899                HistDict[section] = {'sigDict':{},'peaks':[]}
4900                section = 'Index Peak List'
4901                histItems += [section]
4902                HistDict[section] = [[],[]]
4903                section = 'Unit Cells List'
4904                histItems += [section]
4905                HistDict[section] = []
4906                section = 'Reflection Lists'
4907                histItems += [section]
4908                HistDict[section] = {}
4909            # elif 'SASD' in Aname:             
4910            #     section = 'Substances'
4911            #     histItems += [section]
4912            #     HistDict[section] = G2pdG.SetDefaultSubstances()  # this needs to be moved
4913            #     section = 'Sample Parameters'
4914            #     histItems += [section]
4915            #     HistDict[section] = Sample
4916            #     section = 'Models'
4917            #     histItems += [section]
4918            #     HistDict[section] = G2pdG.SetDefaultSASDModel() # this needs to be moved
4919            valuesdict = {
4920                'wtFactor':1.0,'Dummy':False,'ranId':ran.randint(0,sys.maxsize),'Offset':[0.0,0.0],'delOffset':0.02*Ymax,
4921                'refOffset':-0.1*Ymax,'refDelt':0.1*Ymax,'Yminmax':[Ymin,Ymax]}
4922            # if Aname is already in the project replace it
4923            for j in self.proj.names:
4924                if j[0] == Aname: 
4925                    G2fil.G2Print('Replacing "{}" in project'.format(Aname))
4926                    break
4927            else:
4928                G2fil.G2Print('Adding "{}" to project'.format(Aname))
4929                self.proj.names.append([Aname]+
4930                        [u'Comments',u'Limits',u'Background',u'Instrument Parameters',
4931                         u'Sample Parameters', u'Peak List', u'Index Peak List',
4932                         u'Unit Cells List', u'Reflection Lists'])
4933            HistDict['data'] = [valuesdict,
4934                    [np.array(X),np.array(Y),np.array(W),np.zeros(N),np.zeros(N),np.zeros(N)]]
4935            self.proj.data[Aname] = HistDict
4936            IntgOutList.append(self.proj.histogram(Aname))
4937        return IntgOutList
4938
4939##########################
4940# Command Line Interface #
4941##########################
4942# Each of these takes an argparse.Namespace object as their argument,
4943# representing the parsed command-line arguments for the relevant subcommand.
4944# The argument specification for each is in the subcommands dictionary (see
4945# below)
4946
4947commandhelp={}
4948commandhelp["create"] = "creates a GSAS-II project, optionally adding histograms and/or phases"
4949def create(args):
4950    """Implements the create command-line subcommand. This creates a GSAS-II project, optionally adding histograms and/or phases::
4951
4952  usage: GSASIIscriptable.py create [-h] [-d HISTOGRAMS [HISTOGRAMS ...]]
4953                                  [-i IPARAMS [IPARAMS ...]]
4954                                  [-p PHASES [PHASES ...]]
4955                                  filename
4956                                 
4957positional arguments::
4958
4959  filename              the project file to create. should end in .gpx
4960
4961optional arguments::
4962
4963  -h, --help            show this help message and exit
4964  -d HISTOGRAMS [HISTOGRAMS ...], --histograms HISTOGRAMS [HISTOGRAMS ...]
4965                        list of datafiles to add as histograms
4966  -i IPARAMS [IPARAMS ...], --iparams IPARAMS [IPARAMS ...]
4967                        instrument parameter file, must be one for every
4968                        histogram
4969  -p PHASES [PHASES ...], --phases PHASES [PHASES ...]
4970                        list of phases to add. phases are automatically
4971                        associated with all histograms given.
4972
4973    """
4974    proj = G2Project(gpxfile=args.filename)
4975
4976    hist_objs = []
4977    if args.histograms:
4978        for h,i in zip(args.histograms,args.iparams):
4979            G2fil.G2Print("Adding histogram from",h,"with instparm ",i)
4980            hist_objs.append(proj.add_powder_histogram(h, i))
4981
4982    if args.phases: 
4983        for p in args.phases:
4984            G2fil.G2Print("Adding phase from",p)
4985            proj.add_phase(p, histograms=hist_objs)
4986        G2fil.G2Print('Linking phase(s) to histogram(s):')
4987        for h in hist_objs:
4988            G2fil.G2Print ('   '+h.name)
4989
4990    proj.save()
4991
4992commandhelp["add"] = "adds histograms and/or phases to GSAS-II project"
4993def add(args):
4994    """Implements the add command-line subcommand. This adds histograms and/or phases to GSAS-II project::
4995
4996  usage: GSASIIscriptable.py add [-h] [-d HISTOGRAMS [HISTOGRAMS ...]]
4997                               [-i IPARAMS [IPARAMS ...]]
4998                               [-hf HISTOGRAMFORMAT] [-p PHASES [PHASES ...]]
4999                               [-pf PHASEFORMAT] [-l HISTLIST [HISTLIST ...]]
5000                               filename
5001
5002
5003positional arguments::
5004
5005  filename              the project file to open. Should end in .gpx
5006
5007optional arguments::
5008
5009  -h, --help            show this help message and exit
5010  -d HISTOGRAMS [HISTOGRAMS ...], --histograms HISTOGRAMS [HISTOGRAMS ...]
5011                        list of datafiles to add as histograms
5012  -i IPARAMS [IPARAMS ...], --iparams IPARAMS [IPARAMS ...]
5013                        instrument parameter file, must be one for every
5014                        histogram
5015  -hf HISTOGRAMFORMAT, --histogramformat HISTOGRAMFORMAT
5016                        format hint for histogram import. Applies to all
5017                        histograms
5018  -p PHASES [PHASES ...], --phases PHASES [PHASES ...]
5019                        list of phases to add. phases are automatically
5020                        associated with all histograms given.
5021  -pf PHASEFORMAT, --phaseformat PHASEFORMAT
5022                        format hint for phase import. Applies to all phases.
5023                        Example: -pf CIF
5024  -l HISTLIST [HISTLIST ...], --histlist HISTLIST [HISTLIST ...]
5025                        list of histgram indices to associate with added
5026                        phases. If not specified, phases are associated with
5027                        all previously loaded histograms. Example: -l 2 3 4
5028   
5029    """
5030    proj = G2Project(args.filename)
5031
5032    if args.histograms:
5033        for h,i in zip(args.histograms,args.iparams):
5034            G2fil.G2Print("Adding histogram from",h,"with instparm ",i)
5035            proj.add_powder_histogram(h, i, fmthint=args.histogramformat)
5036
5037    if args.phases: 
5038        if not args.histlist:
5039            histlist = proj.histograms()
5040        else:
5041            histlist = [proj.histogram(i) for i in args.histlist]
5042
5043        for p in args.phases:
5044            G2fil.G2Print("Adding phase from",p)
5045            proj.add_phase(p, histograms=histlist, fmthint=args.phaseformat)
5046           
5047        if not args.histlist:
5048            G2fil.G2Print('Linking phase(s) to all histogram(s)')
5049        else:
5050            G2fil.G2Print('Linking phase(s) to histogram(s):')
5051            for h in histlist:
5052                G2fil.G2Print('   '+h.name)
5053    proj.save()
5054
5055
5056commandhelp["dump"] = "Shows the contents of a GSAS-II project"
5057def dump(args):
5058    """Implements the dump command-line subcommand, which shows the contents of a GSAS-II project::
5059
5060       usage: GSASIIscriptable.py dump [-h] [-d] [-p] [-r] files [files ...]
5061
5062positional arguments::
5063
5064  files
5065
5066optional arguments::
5067
5068  -h, --help        show this help message and exit
5069  -d, --histograms  list histograms in files, overrides --raw
5070  -p, --phases      list phases in files, overrides --raw
5071  -r, --raw         dump raw file contents, default
5072 
5073    """
5074    if not args.histograms and not args.phases:
5075        args.raw = True
5076    if args.raw:
5077        import IPython.lib.pretty as pretty
5078
5079    for fname in args.files:
5080        if args.raw:
5081            proj, nameList = LoadDictFromProjFile(fname)
5082            print("file:", fname)
5083            print("names:", nameList)
5084            for key, val in proj.items():
5085                print(key, ":")
5086                pretty.pprint(val)
5087        else:
5088            proj = G2Project(fname)
5089            if args.histograms:
5090                hists = proj.histograms()
5091                for h in hists:
5092                    print(fname, "hist", h.id, h.name)
5093            if args.phases:
5094                phase_list = proj.phases()
5095                for p in phase_list:
5096                    print(fname, "phase", p.id, p.name)
5097
5098
5099commandhelp["browse"] = "Load a GSAS-II project and then open a IPython shell to browse it"
5100def IPyBrowse(args):
5101    """Load a .gpx file and then open a IPython shell to browse it::
5102
5103  usage: GSASIIscriptable.py browse [-h] files [files ...]
5104
5105positional arguments::
5106
5107  files       list of files to browse
5108
5109optional arguments::
5110
5111  -h, --help  show this help message and exit
5112
5113    """
5114    for fname in args.files:
5115        proj, nameList = LoadDictFromProjFile(fname)
5116        msg = "\nfname {} loaded into proj (dict) with names in nameList".format(fname)
5117        GSASIIpath.IPyBreak_base(msg)
5118        break
5119
5120
5121commandhelp["refine"] = '''
5122Conducts refinements on GSAS-II projects according to a list of refinement
5123steps in a JSON dict
5124'''
5125def refine(args):
5126    """Implements the refine command-line subcommand:
5127    conducts refinements on GSAS-II projects according to a JSON refinement dict::
5128
5129        usage: GSASIIscriptable.py refine [-h] gpxfile [refinements]
5130
5131positional arguments::
5132
5133  gpxfile      the project file to refine
5134  refinements  json file of refinements to apply. if not present refines file
5135               as-is
5136
5137optional arguments::
5138
5139  -h, --help   show this help message and exit
5140 
5141    """
5142    proj = G2Project(args.gpxfile)
5143    if args.refinements is None:
5144        proj.refine()
5145    else:
5146        import json
5147        with open(args.refinements) as refs:
5148            refs = json.load(refs)
5149        if type(refs) is not dict:
5150            raise G2ScriptException("Error: JSON object must be a dict.")
5151        if "code" in refs:
5152            print("executing code:\n|  ",'\n|  '.join(refs['code']))
5153            exec('\n'.join(refs['code']))
5154        proj.do_refinements(refs['refinements'])
5155
5156
5157commandhelp["export"] = "Export phase as CIF"
5158def export(args):
5159    """Implements the export command-line subcommand: Exports phase as CIF::
5160
5161      usage: GSASIIscriptable.py export [-h] gpxfile phase exportfile
5162
5163positional arguments::
5164
5165  gpxfile     the project file from which to export
5166  phase       identifier of phase to export
5167  exportfile  the .cif file to export to
5168
5169optional arguments::
5170
5171  -h, --help  show this help message and exit
5172
5173    """
5174    proj = G2Project(args.gpxfile)
5175    phase = proj.phase(args.phase)
5176    phase.export_CIF(args.exportfile)
5177
5178
5179def _args_kwargs(*args, **kwargs):
5180    return args, kwargs
5181
5182# A dictionary of the name of each subcommand, and a tuple
5183# of its associated function and a list of its arguments
5184# The arguments are passed directly to the add_argument() method
5185# of an argparse.ArgumentParser
5186
5187subcommands = {"create":
5188               (create, [_args_kwargs('filename',
5189                                      help='the project file to create. should end in .gpx'),
5190
5191                         _args_kwargs('-d', '--histograms',
5192                                      nargs='+',
5193                                      help='list of datafiles to add as histograms'),
5194                                     
5195                         _args_kwargs('-i', '--iparams',
5196                                      nargs='+',
5197                                      help='instrument parameter file, must be one'
5198                                           ' for every histogram'
5199                                      ),
5200
5201                         _args_kwargs('-p', '--phases',
5202                                      nargs='+',
5203                                      help='list of phases to add. phases are '
5204                                           'automatically associated with all '
5205                                           'histograms given.')]),
5206               "add": (add, [_args_kwargs('filename',
5207                                      help='the project file to open. Should end in .gpx'),
5208
5209                         _args_kwargs('-d', '--histograms',
5210                                      nargs='+',
5211                                      help='list of datafiles to add as histograms'),
5212                                     
5213                         _args_kwargs('-i', '--iparams',
5214                                      nargs='+',
5215                                      help='instrument parameter file, must be one'
5216                                           ' for every histogram'
5217                                      ),
5218                                     
5219                         _args_kwargs('-hf', '--histogramformat',
5220                                      help='format hint for histogram import. Applies to all'
5221                                           ' histograms'
5222                                      ),
5223
5224                         _args_kwargs('-p', '--phases',
5225                                      nargs='+',
5226                                      help='list of phases to add. phases are '
5227                                           'automatically associated with all '
5228                                           'histograms given.'),
5229
5230                         _args_kwargs('-pf', '--phaseformat',
5231                                      help='format hint for phase import. Applies to all'
5232                                           ' phases. Example: -pf CIF'
5233                                      ),
5234                                     
5235                         _args_kwargs('-l', '--histlist',
5236                                      nargs='+',
5237                                      help='list of histgram indices to associate with added'
5238                                           ' phases. If not specified, phases are'
5239                                           ' associated with all previously loaded'
5240                                           ' histograms. Example: -l 2 3 4')]),
5241                                           
5242               "dump": (dump, [_args_kwargs('-d', '--histograms',
5243                                     action='store_true',
5244                                     help='list histograms in files, overrides --raw'),
5245
5246                               _args_kwargs('-p', '--phases',
5247                                            action='store_true',
5248                                            help='list phases in files, overrides --raw'),
5249
5250                               _args_kwargs('-r', '--raw',
5251                                      action='store_true', help='dump raw file contents, default'),
5252
5253                               _args_kwargs('files', nargs='+')]),
5254
5255               "refine":
5256               (refine, [_args_kwargs('gpxfile', help='the project file to refine'),
5257                         _args_kwargs('refinements',
5258                                      help='JSON file of refinements to apply. if not present'
5259                                           ' refines file as-is',
5260                                      default=None,
5261                                      nargs='?')]),
5262
5263               "export": (export, [_args_kwargs('gpxfile',
5264                                                help='the project file from which to export'),
5265                                   _args_kwargs('phase', help='identifier of phase to export'),
5266                                   _args_kwargs('exportfile', help='the .cif file to export to')]),
5267               "browse": (IPyBrowse, [_args_kwargs('files', nargs='+',
5268                                                   help='list of files to browse')])}
5269
5270
5271def main():
5272    '''The command-line interface for calling GSASIIscriptable as a shell command,
5273    where it is expected to be called as::
5274
5275       python GSASIIscriptable.py <subcommand> <file.gpx> <options>
5276
5277    The following subcommands are defined:
5278
5279        * create, see :func:`create`
5280        * add, see :func:`add`
5281        * dump, see :func:`dump`
5282        * refine, see :func:`refine`
5283        * export, :func:`export`
5284        * browse, see :func:`IPyBrowse`
5285
5286    .. seealso::
5287        :func:`create`
5288        :func:`add`
5289        :func:`dump`
5290        :func:`refine`
5291        :func:`export`
5292        :func:`IPyBrowse`
5293    '''
5294    parser = argparse.ArgumentParser(description=
5295        "Use of "+os.path.split(__file__)[1]+" Allows GSAS-II actions from command line."
5296        )
5297    subs = parser.add_subparsers()
5298
5299    # Create all of the specified subparsers
5300    for name, (func, args) in subcommands.items():
5301        new_parser = subs.add_parser(name,help=commandhelp.get(name),
5302                                     description='Command "'+name+'" '+commandhelp.get(name))
5303        for listargs, kwds in args:
5304            new_parser.add_argument(*listargs, **kwds)
5305        new_parser.set_defaults(func=func)
5306
5307    # Parse and trigger subcommand
5308    result = parser.parse_args()
5309    result.func(result)
5310
5311if __name__ == '__main__':
5312    #fname='/tmp/corundum-template.gpx'
5313    #prj = G2Project(fname)
5314    main()
Note: See TracBrowser for help on using the repository browser.