source: trunk/GSASIIscriptable.py @ 4031

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

introduce external autoInt

  • 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-18 23:21:39 +0000 (Tue, 18 Jun 2019) $
5# $Author: toby $
6# $Revision: 4031 $
7# $URL: trunk/GSASIIscriptable.py $
8# $Id: GSASIIscriptable.py 4031 2019-06-18 23:21:39Z 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 PDFobj: A :class:`G2PDF` object which may be
2109          in a separate project or the dict associated with the
2110          PDF object (G2PDF.data).
2111        :param histogram: A reference to a histogram,
2112          which can be reference by object, name, or number.
2113        :returns: A :class:`G2PDF` object for the PDF entry
2114        '''
2115        LoadG2fil()
2116        PDFname = 'PDF ' + self.histogram(histogram).name[4:]
2117        PDFdict = {'data':None}
2118        for i in 'PDF Controls', 'PDF Peaks':
2119            PDFdict[i] = copy.deepcopy(PDFobj[i])
2120        self.names.append([PDFname]+['PDF Controls', 'PDF Peaks'])
2121        self.data[PDFname] = PDFdict
2122        for i in 'I(Q)','S(Q)','F(Q)','G(R)':
2123            self.data[PDFname]['PDF Controls'][i] = []
2124        G2fil.G2Print('Adding "{}" to project'.format(PDFname))
2125        return G2PDF(self.data[PDFname], PDFname, self)
2126       
2127    def add_PDF(self, prmfile, histogram):
2128        '''Creates a PDF entry that can be used to compute a PDF.
2129        Note that this command places an entry in the project,
2130        but :meth:`G2PDF.calculate` must be used to actually perform
2131        the computation.
2132
2133        :param str datafile: The powder data file to read, a filename.
2134        :param histogram: A reference to a histogram,
2135          which can be reference by object, name, or number.
2136        :returns: A :class:`G2PDF` object for the PDF entry
2137        '''
2138       
2139        LoadG2fil()
2140        PDFname = 'PDF ' + self.histogram(histogram).name[4:]
2141        peaks = {'Limits':[1.,5.],'Background':[2,[0.,-0.2*np.pi],False],'Peaks':[]}
2142        Controls = {
2143            'Sample':{'Name':self.histogram(histogram).name,'Mult':1.0},
2144            'Sample Bkg.':{'Name':'','Mult':-1.0,'Refine':False},
2145            'Container':{'Name':'','Mult':-1.0,'Refine':False},
2146            'Container Bkg.':{'Name':'','Mult':-1.0},
2147            'ElList':{},
2148            'Geometry':'Cylinder','Diam':1.0,'Pack':0.50,'Form Vol':0.0,'Flat Bkg':0,
2149            'DetType':'Area detector','ObliqCoeff':0.2,'Ruland':0.025,'QScaleLim':[20,25],
2150            'Lorch':False,'BackRatio':0.0,'Rmax':100.,'noRing':False,'IofQmin':1.0,'Rmin':1.0,
2151            'I(Q)':[],'S(Q)':[],'F(Q)':[],'G(R)':[]}
2152
2153        fo = open(prmfile,'r')
2154        S = fo.readline()
2155        while S:
2156            if '#' in S:
2157                S = fo.readline()
2158                continue
2159            key,val = S.split(':',1)
2160            try:
2161                Controls[key] = eval(val)
2162            except:
2163                Controls[key] = val.strip()
2164            S = fo.readline()
2165        fo.close()
2166        Controls['Sample']['Name'] = self.histogram(histogram).name
2167        for i in 'Sample Bkg.','Container','Container Bkg.':
2168            Controls[i]['Name'] = ''
2169        PDFdict = {'data':None,'PDF Controls':Controls, 'PDF Peaks':peaks}
2170        self.names.append([PDFname]+['PDF Controls', 'PDF Peaks'])
2171        self.data[PDFname] = PDFdict
2172        G2fil.G2Print('Adding "{}" to project'.format(PDFname))
2173        return G2PDF(self.data[PDFname], PDFname, self)
2174
2175    def seqref(self):
2176        """
2177        Returns a sequential refinement results object, if present
2178
2179        :returns: A :class:`G2SeqRefRes` object or None if not present
2180        """
2181        if 'Sequential results' not in self.data: return
2182        return G2SeqRefRes(self.data['Sequential results']['data'], self)
2183   
2184    def update_ids(self):
2185        """Makes sure all phases and histograms have proper hId and pId"""
2186        # Translated from GetUsedHistogramsAndPhasesfromTree,
2187        #   GSASIIdataGUI.py:4107
2188        for i, h in enumerate(self.histograms()):
2189            h.id = i
2190        for i, p in enumerate(self.phases()):
2191            p.id = i
2192
2193    def do_refinements(self, refinements, histogram='all', phase='all',
2194                       outputnames=None, makeBack=False):
2195        """Conducts one or a series of refinements according to the
2196           input provided in parameter refinements. This is a wrapper
2197           around :meth:`iter_refinements`
2198
2199        :param list refinements: A list of dictionaries specifiying changes to be made to
2200            parameters before refinements are conducted.
2201            See the :ref:`Refinement_recipe` section for how this is defined.
2202        :param str histogram: Name of histogram for refinements to be applied
2203            to, or 'all'; note that this can be overridden for each refinement
2204            step via a "histograms" entry in the dict.
2205        :param str phase: Name of phase for refinements to be applied to, or
2206            'all'; note that this can be overridden for each refinement
2207            step via a "phases" entry in the dict.
2208        :param list outputnames: Provides a list of project (.gpx) file names
2209            to use for each refinement step (specifying None skips the save step).
2210            See :meth:`save`.
2211            Note that this can be overridden using an "output" entry in the dict.
2212        :param bool makeBack: determines if a backup ).bckX.gpx) file is made
2213            before a refinement is performed. The default is False.
2214           
2215        To perform a single refinement without changing any parameters, use this
2216        call:
2217
2218        .. code-block::  python
2219       
2220            my_project.do_refinements([])
2221        """
2222       
2223        for proj in self.iter_refinements(refinements, histogram, phase,
2224                                          outputnames, makeBack):
2225            pass
2226        return self
2227
2228    def iter_refinements(self, refinements, histogram='all', phase='all',
2229                         outputnames=None, makeBack=False):
2230        """Conducts a series of refinements, iteratively. Stops after every
2231        refinement and yields this project, to allow error checking or
2232        logging of intermediate results. Parameter use is the same as for
2233        :meth:`do_refinements` (which calls this method).
2234
2235        >>> def checked_refinements(proj):
2236        ...     for p in proj.iter_refinements(refs):
2237        ...         # Track intermediate results
2238        ...         log(p.histogram('0').residuals)
2239        ...         log(p.phase('0').get_cell())
2240        ...         # Check if parameter diverged, nonsense answer, or whatever
2241        ...         if is_something_wrong(p):
2242        ...             raise Exception("I need a human!")
2243
2244           
2245        """
2246        if outputnames:
2247            if len(refinements) != len(outputnames):
2248                raise ValueError("Should have same number of outputs to"
2249                                 "refinements")
2250        else:
2251            outputnames = [None for r in refinements]
2252
2253        for output, refinedict in zip(outputnames, refinements):
2254            if 'histograms' in refinedict:
2255                hist = refinedict['histograms']
2256            else:
2257                hist = histogram
2258            if 'phases' in refinedict:
2259                ph = refinedict['phases']
2260            else:
2261                ph = phase
2262            if 'output' in refinedict:
2263                output = refinedict['output']
2264            self.set_refinement(refinedict, hist, ph)
2265            # Handle 'once' args - refinements that are disabled after this
2266            # refinement
2267            if 'once' in refinedict:
2268                temp = {'set': refinedict['once']}
2269                self.set_refinement(temp, hist, ph)
2270
2271            if output:
2272                self.save(output)
2273
2274            if 'skip' not in refinedict:
2275                self.refine(makeBack=makeBack)
2276            yield self
2277
2278            # Handle 'once' args - refinements that are disabled after this
2279            # refinement
2280            if 'once' in refinedict:
2281                temp = {'clear': refinedict['once']}
2282                self.set_refinement(temp, hist, ph)
2283            if 'call' in refinedict:
2284                fxn = refinedict['call']
2285                if callable(fxn):
2286                    fxn(*refinedict.get('callargs',[self]))
2287                elif callable(eval(fxn)):
2288                    eval(fxn)(*refinedict.get('callargs',[self]))
2289                else:
2290                    raise G2ScriptException("Error: call value {} is not callable".format(fxn))
2291
2292    def set_refinement(self, refinement, histogram='all', phase='all'):
2293        """Apply specified refinements to a given histogram(s) or phase(s).
2294
2295        :param dict refinement: The refinements to be conducted
2296        :param histogram: Specifies either 'all' (default), a single histogram or
2297          a list of histograms. Histograms may be specified as histogram objects
2298          (see :class:`G2PwdrData`), the histogram name (str) or the index number (int)
2299          of the histogram in the project, numbered starting from 0.
2300          Omitting the parameter or the string 'all' indicates that parameters in
2301          all histograms should be set.
2302        :param phase: Specifies either 'all' (default), a single phase or
2303          a list of phases. Phases may be specified as phase objects
2304          (see :class:`G2Phase`), the phase name (str) or the index number (int)
2305          of the phase in the project, numbered starting from 0.
2306          Omitting the parameter or the string 'all' indicates that parameters in
2307          all phases should be set.
2308
2309        Note that refinement parameters are categorized as one of three types:
2310
2311        1. Histogram parameters
2312        2. Phase parameters
2313        3. Histogram-and-Phase (HAP) parameters
2314       
2315        .. seealso::
2316            :meth:`G2PwdrData.set_refinements`
2317            :meth:`G2PwdrData.clear_refinements`
2318            :meth:`G2Phase.set_refinements`
2319            :meth:`G2Phase.clear_refinements`
2320            :meth:`G2Phase.set_HAP_refinements`
2321            :meth:`G2Phase.clear_HAP_refinements`"""
2322
2323        if histogram == 'all':
2324            hists = self.histograms()
2325        elif isinstance(histogram, list) or isinstance(histogram, tuple):
2326            hists = []
2327            for h in histogram:
2328                if isinstance(h, str) or isinstance(h, int):
2329                    hists.append(self.histogram(h))
2330                else:
2331                    hists.append(h)
2332        elif isinstance(histogram, str) or isinstance(histogram, int):
2333            hists = [self.histogram(histogram)]
2334        else:
2335            hists = [histogram]
2336
2337        if phase == 'all':
2338            phases = self.phases()
2339        elif isinstance(phase, list) or isinstance(phase, tuple):
2340            phases = []
2341            for ph in phase:
2342                if isinstance(ph, str) or isinstance(ph, int):
2343                    phases.append(self.phase(ph))
2344                else:
2345                    phases.append(ph)
2346        elif isinstance(phase, str) or isinstance(phase, int):
2347            phases = [self.phase(phase)]
2348        else:
2349            phases = [phase]
2350
2351        pwdr_set = {}
2352        phase_set = {}
2353        hap_set = {}
2354        for key, val in refinement.get('set', {}).items():
2355            # Apply refinement options
2356            if G2PwdrData.is_valid_refinement_key(key):
2357                pwdr_set[key] = val
2358            elif G2Phase.is_valid_refinement_key(key):
2359                phase_set[key] = val
2360            elif G2Phase.is_valid_HAP_refinement_key(key):
2361                hap_set[key] = val
2362            else:
2363                raise ValueError("Unknown refinement key", key)
2364
2365        for hist in hists:
2366            hist.set_refinements(pwdr_set)
2367        for phase in phases:
2368            phase.set_refinements(phase_set)
2369        for phase in phases:
2370            phase.set_HAP_refinements(hap_set, hists)
2371
2372        pwdr_clear = {}
2373        phase_clear = {}
2374        hap_clear = {}
2375        for key, val in refinement.get('clear', {}).items():
2376            # Clear refinement options
2377            if G2PwdrData.is_valid_refinement_key(key):
2378                pwdr_clear[key] = val
2379            elif G2Phase.is_valid_refinement_key(key):
2380                phase_clear[key] = val
2381            elif G2Phase.is_valid_HAP_refinement_key(key):
2382                hap_set[key] = val
2383            else:
2384                raise ValueError("Unknown refinement key", key)
2385
2386        for hist in hists:
2387            hist.clear_refinements(pwdr_clear)
2388        for phase in phases:
2389            phase.clear_refinements(phase_clear)
2390        for phase in phases:
2391            phase.clear_HAP_refinements(hap_clear, hists)
2392
2393    def index_ids(self):
2394        self.save()
2395        return G2strIO.GetUsedHistogramsAndPhases(self.filename)
2396
2397    def add_constraint_raw(self, cons_scope, constr):
2398        """Adds a constraint of type consType to the project.
2399        cons_scope should be one of "Hist", "Phase", "HAP", or "Global".
2400
2401        WARNING it does not check the constraint is well-constructed"""
2402        constrs = self.data['Constraints']['data']
2403        if 'Global' not in constrs:
2404            constrs['Global'] = []
2405        constrs[cons_scope].append(constr)
2406
2407    def hold_many(self, vars, type):
2408        """Apply holds for all the variables in vars, for constraint of a given type.
2409
2410        type is passed directly to add_constraint_raw as consType
2411
2412        :param list vars: A list of variables to hold. Either :class:`GSASIIobj.G2VarObj` objects,
2413            string variable specifiers, or arguments for :meth:`make_var_obj`
2414        :param str type: A string constraint type specifier. See
2415            :class:`G2Project.add_constraint_raw`
2416
2417        """
2418        for var in vars:
2419            if isinstance(var, str):
2420                var = self.make_var_obj(var)
2421            elif not isinstance(var, G2obj.G2VarObj):
2422                var = self.make_var_obj(*var)
2423            self.add_constraint_raw(type, [[1.0, var], None, None, 'h'])
2424
2425    def make_var_obj(self, phase=None, hist=None, varname=None, atomId=None,
2426                     reloadIdx=True):
2427        """Wrapper to create a G2VarObj. Takes either a string representation ("p:h:name:a")
2428        or individual names of phase, histogram, varname, and atomId.
2429
2430        Automatically converts string phase, hist, or atom names into the ID required
2431        by G2VarObj."""
2432
2433        if reloadIdx:
2434            self.index_ids()
2435
2436        # If string representation, short circuit
2437        if hist is None and varname is None and atomId is None:
2438            if isinstance(phase, str) and ':' in phase:
2439                return G2obj.G2VarObj(phase)
2440
2441        # Get phase index
2442        phaseObj = None
2443        if isinstance(phase, G2Phase):
2444            phaseObj = phase
2445            phase = G2obj.PhaseRanIdLookup[phase.ranId]
2446        elif phase in self.data['Phases']:
2447            phaseObj = self.phase(phase)
2448            phaseRanId = phaseObj.ranId
2449            phase = G2obj.PhaseRanIdLookup[phaseRanId]
2450
2451        # Get histogram index
2452        if isinstance(hist, G2PwdrData):
2453            hist = G2obj.HistRanIdLookup[hist.ranId]
2454        elif hist in self.data:
2455            histRanId = self.histogram(hist).ranId
2456            hist = G2obj.HistRanIdLookup[histRanId]
2457
2458        # Get atom index (if any)
2459        if isinstance(atomId, G2AtomRecord):
2460            atomId = G2obj.AtomRanIdLookup[phase][atomId.ranId]
2461        elif phaseObj:
2462            atomObj = phaseObj.atom(atomId)
2463            if atomObj:
2464                atomRanId = atomObj.ranId
2465                atomId = G2obj.AtomRanIdLookup[phase][atomRanId]
2466
2467        return G2obj.G2VarObj(phase, hist, varname, atomId)
2468
2469    def add_image(self, imagefile, fmthint=None, defaultImage=None):
2470        """Load an image into a project
2471
2472        :param str imagefile: The image file to read, a filename.
2473        :param str fmthint: If specified, only importers where the format name
2474          (reader.formatName, as shown in Import menu) contains the
2475          supplied string will be tried as importers. If not specified, all
2476          importers consistent with the file extension will be tried
2477          (equivalent to "guess format" in menu).
2478        :param str defaultImage: The name of an image to use as a default for
2479          setting parameters for the image file to read.
2480
2481        :returns: a list of G2Image object for the added image(s) [this routine
2482          has not yet been tested with files with multiple images...]
2483        """
2484        LoadG2fil()
2485        imagefile = os.path.abspath(os.path.expanduser(imagefile))
2486        readers = import_generic(imagefile, Readers['Image'], fmthint=fmthint)
2487        objlist = []
2488        for rd in readers:
2489            if rd.SciPy:        #was default read by scipy; needs 1 time fixes
2490                G2fil.G2Print('Warning: Image {} read by generic SciPy import. Image parameters likely wrong'.format(imagefile))
2491                #see this: G2IO.EditImageParms(self,rd.Data,rd.Comments,rd.Image,imagefile)
2492                rd.SciPy = False
2493            rd.readfilename = imagefile
2494            if rd.repeatcount == 1 and not rd.repeat: # skip image number if only one in set
2495                rd.Data['ImageTag'] = None
2496            else:
2497                rd.Data['ImageTag'] = rd.repeatcount
2498            rd.Data['formatName'] = rd.formatName
2499            if rd.sumfile:
2500                rd.readfilename = rd.sumfile
2501            # Load generic metadata, as configured
2502            G2fil.GetColumnMetadata(rd)
2503            # Code from G2IO.LoadImage2Tree(rd.readfilename,self,rd.Comments,rd.Data,rd.Npix,rd.Image)
2504            Imax = np.amax(rd.Image)
2505            ImgNames = [i[0] for i in self.names if i[0].startswith('IMG ')]
2506            TreeLbl = 'IMG '+os.path.basename(imagefile)
2507            ImageTag = rd.Data.get('ImageTag')
2508            if ImageTag:
2509                TreeLbl += ' #'+'%04d'%(ImageTag)
2510                imageInfo = (imagefile,ImageTag)
2511            else:
2512                imageInfo = imagefile
2513            TreeName = G2obj.MakeUniqueLabel(TreeLbl,ImgNames)
2514            # MT dict to contain image info
2515            ImgDict = {}
2516            ImgDict['data'] = [rd.Npix,imageInfo]
2517            ImgDict['Comments'] = rd.Comments
2518            if defaultImage:
2519                if isinstance(defaultImage, G2Image):
2520                    if defaultImage.proj == self:
2521                        defaultImage = self.data[defaultImage.name]['data']
2522                    else:
2523                        raise Exception("Image {} not in current selected project".format(defaultImage.name))
2524                elif defaultImage in self._images():
2525                    defaultImage = self.data[defaultImage]['data']
2526                else:
2527                    defaultImage = None
2528            Data = rd.Data
2529            if defaultImage:
2530                Data = copy.copy(defaultImage)
2531                Data['showLines'] = True
2532                Data['ring'] = []
2533                Data['rings'] = []
2534                Data['cutoff'] = 10.
2535                Data['pixLimit'] = 20
2536                Data['edgemin'] = 100000000
2537                Data['calibdmin'] = 0.5
2538                Data['calibskip'] = 0
2539                Data['ellipses'] = []
2540                Data['calibrant'] = ''
2541                Data['GonioAngles'] = [0.,0.,0.]
2542                Data['DetDepthRef'] = False
2543            else:
2544                Data['type'] = 'PWDR'
2545                Data['color'] = GSASIIpath.GetConfigValue('Contour_color','Paired')
2546                if 'tilt' not in Data:          #defaults if not preset in e.g. Bruker importer
2547                    Data['tilt'] = 0.0
2548                    Data['rotation'] = 0.0
2549                    Data['pixLimit'] = 20
2550                    Data['calibdmin'] = 0.5
2551                    Data['cutoff'] = 10.
2552                Data['showLines'] = False
2553                Data['calibskip'] = 0
2554                Data['ring'] = []
2555                Data['rings'] = []
2556                Data['edgemin'] = 100000000
2557                Data['ellipses'] = []
2558                Data['GonioAngles'] = [0.,0.,0.]
2559                Data['DetDepth'] = 0.
2560                Data['DetDepthRef'] = False
2561                Data['calibrant'] = ''
2562                Data['IOtth'] = [5.0,50.0]
2563                Data['LRazimuth'] = [0.,180.]
2564                Data['azmthOff'] = 0.0
2565                Data['outChannels'] = 2500
2566                Data['outAzimuths'] = 1
2567                Data['centerAzm'] = False
2568                Data['fullIntegrate'] = GSASIIpath.GetConfigValue('fullIntegrate',True)
2569                Data['setRings'] = False
2570                Data['background image'] = ['',-1.0]                           
2571                Data['dark image'] = ['',-1.0]
2572                Data['Flat Bkg'] = 0.0
2573                Data['Oblique'] = [0.5,False]
2574            Data['setDefault'] = False
2575            Data['range'] = [(0,Imax),[0,Imax]]
2576            ImgDict['Image Controls'] = Data
2577            ImgDict['Masks'] = {'Points':[],'Rings':[],'Arcs':[],'Polygons':[],
2578                                'Frames':[],'Thresholds':[(0,Imax),[0,Imax]]}
2579            ImgDict['Stress/Strain']  = {'Type':'True','d-zero':[],'Sample phi':0.0,
2580                                             'Sample z':0.0,'Sample load':0.0}
2581            self.names.append([TreeName]+['Comments','Image Controls','Masks','Stress/Strain'])
2582            self.data[TreeName] = ImgDict
2583            del rd.Image
2584            objlist.append(G2Image(self.data[TreeName], TreeName, self))
2585        return objlist
2586   
2587    def get_Controls(self, control):
2588        '''Return project controls settings
2589
2590        :param str control: the item to be returned. See below for allowed values.
2591        :returns: The value for the control.
2592
2593        Allowed values for parameter control:
2594
2595            * cycles: the maximum number of cycles (returns int)
2596            * sequential: the histograms used for a sequential refinement as a list
2597              of histogram names or an empty list when in non-sequential mode.
2598            * seqCopy: returns True or False. True indicates that results from
2599              each sequential fit are used as the starting point for the next
2600              histogram.
2601            * Anything else returns the value in the Controls dict, if present. An
2602              exception is raised if the control value is not present.
2603
2604        .. seealso::
2605            :meth:`set_Controls`
2606        '''
2607        if control == 'cycles':
2608            return self.data['Controls']['data']['max cyc']
2609        elif control == 'sequential':
2610            return self.data['Controls']['data']['Seq Data']
2611        elif control in self.data['Controls']['data']:
2612            return self.data['Controls']['data'][control]
2613        elif control == 'seqCopy':
2614            return self.data['Controls']['data']['Copy2Next']
2615        else:
2616            G2fil.G2Print('Defined Controls:',self.data['Controls']['data'].keys())
2617            raise Exception('{} is an invalid control value'.format(control))
2618       
2619    def set_Controls(self, control, value):
2620        '''Set project controls
2621
2622        :param str control: the item to be set. See below for allowed values.
2623        :param value: the value to be set.
2624
2625        Allowed values for parameter control:
2626
2627            * cycles: sets the maximum number of cycles (value must be int)
2628            * sequential: sets the histograms to be used for a sequential refinement.
2629              Use an empty list to turn off sequential fitting.
2630              The values in the list may be the name of the histogram (a str), or
2631              a ranId or index (int values), see :meth:`histogram`.
2632            * seqCopy: when True, the results from each sequential fit are used as
2633              the starting point for the next. After each fit is is set to False.
2634              Ignored for non-sequential fits.
2635
2636        .. seealso::
2637            :meth:`get_Controls`
2638        '''
2639        if control == 'cycles':
2640            self.data['Controls']['data']['max cyc'] = int(value)
2641        elif control == 'seqCopy':
2642            self.data['Controls']['data']['Copy2Next'] = bool(value)
2643        elif control == 'sequential':
2644            histlist = []
2645            for i,j in enumerate(value):
2646                h = self.histogram(j)
2647                if h:
2648                    histlist.append(h.name)
2649                else:
2650                    raise Exception('item #{} ({}) is an invalid histogram value'
2651                                        .format(i,j))
2652            self.data['Controls']['data']['Seq Data'] = histlist
2653        else:
2654            raise Exception('{} is an invalid control value'.format(control))
2655       
2656    def copyHistParms(self,sourcehist,targethistlist='all',modelist='all'):
2657        '''Copy histogram information from one histogram to others
2658
2659        :param sourcehist: is a histogram object (:class:`G2PwdrData`) or
2660            a histogram name or the index number of the histogram
2661        :param list targethistlist: a list of histograms where each item in the
2662            list can be a histogram object (:class:`G2PwdrData`),
2663            a histogram name or the index number of the histogram.
2664            if the string 'all' (default value), then all histograms in
2665            the project are used.
2666        :param list modelist: May be a list of sections to copy, which may
2667           include 'Background', 'Instrument Parameters', 'Limits' and
2668           'Sample Parameters' (items may be shortened to uniqueness and
2669           capitalization is ignored, so ['b','i','L','s'] will work.)
2670           The default value, 'all' causes the listed sections to
2671
2672        '''
2673        sections = ('Background','Instrument Parameters','Limits',
2674                        'Sample Parameters')
2675        hist_in = self.histogram(sourcehist)
2676        if not hist_in:
2677            raise Exception('{} is not a valid histogram'.format(sourcehist))
2678        if targethistlist == "all":
2679            targethistlist = self.histograms()
2680        if 'all' in modelist:
2681            copysections = sections
2682        else:
2683            copysections = set()
2684            for s in sections:
2685                for m in modelist:
2686                    if s.lower().startswith(m.lower()):
2687                        copysections.add(s)
2688        for h in targethistlist:
2689            hist_out = self.histogram(h)
2690            if not hist_out:
2691                raise Exception('{} is not a valid histogram'.format(h))
2692            for key in copysections: 
2693                hist_out[key] = copy.deepcopy(hist_in[key])
2694
2695    def get_VaryList(self):
2696        '''Returns a list of the refined variables in the
2697        last refinement cycle
2698
2699        :returns: a list of variables or None if no refinement has been
2700          performed.
2701        '''
2702        try:
2703            return self['Covariance']['data']['varyList']
2704        except:
2705            return
2706
2707    def get_ParmList(self):
2708        '''Returns a list of all the parameters defined in the
2709        last refinement cycle
2710
2711        :returns: a list of parameters or None if no refinement has been
2712          performed.
2713        '''
2714        try:
2715            return list(self['Covariance']['data']['parmDict'].keys())
2716        except:
2717            return
2718       
2719    def get_Variable(self,var):
2720        '''Returns the value and standard uncertainty (esd) for a variable
2721        parameters, as defined in the last refinement cycle
2722
2723        :param str var: a variable name of form '<p>:<h>:<name>', such as
2724          ':0:Scale'
2725        :returns: (value,esd) if the parameter is refined or
2726          (value, None) if the variable is in a constraint or is not
2727          refined or None if the parameter is not found.
2728        '''
2729        if var not in self['Covariance']['data']['parmDict']:
2730            return None
2731        val = self['Covariance']['data']['parmDict'][var]
2732        try:
2733            pos = self['Covariance']['data']['varyList'].index(var)
2734            esd = np.sqrt(self['Covariance']['data']['covMatrix'][pos,pos])
2735            return (val,esd)
2736        except ValueError:
2737            return (val,None)
2738
2739    def get_Covariance(self,varList):
2740        '''Returns the values and covariance matrix for a series of variable
2741        parameters. as defined in the last refinement cycle
2742
2743        :param tuple varList: a list of variable names of form '<p>:<h>:<name>'
2744        :returns: (valueList,CovMatrix) where valueList contains the (n) values
2745          in the same order as varList (also length n) and CovMatrix is a
2746          (n x n) matrix. If any variable name is not found in the varyList
2747          then None is returned.
2748
2749        Use this code, where sig provides standard uncertainties for
2750        parameters and where covArray provides the correlation between
2751        off-diagonal terms::
2752
2753            sig = np.sqrt(np.diag(covMatrix))
2754            xvar = np.outer(sig,np.ones_like(sig))
2755            covArray = np.divide(np.divide(covMatrix,xvar),xvar.T)
2756
2757        '''
2758        missing = [i for i in varList if i not in self['Covariance']['data']['varyList']]
2759        if missing:
2760            G2fil.G2Print('Warning: Variable(s) {} were not found in the varyList'.format(missing))
2761            return None
2762        vals = [self['Covariance']['data']['parmDict'][i] for i in varList]
2763        import GSASIImath as G2mth
2764        cov = G2mth.getVCov(varList,
2765                            self['Covariance']['data']['varyList'],
2766                            self['Covariance']['data']['covMatrix'])
2767        return (vals,cov)
2768
2769class G2AtomRecord(G2ObjectWrapper):
2770    """Wrapper for an atom record. Has convenient accessors via @property.
2771
2772    Available accessors: label, type, refinement_flags, coordinates,
2773        occupancy, ranId, id, adp_flag, uiso
2774
2775    Example:
2776
2777    >>> atom = some_phase.atom("O3")
2778    >>> # We can access the underlying data format
2779    >>> atom.data
2780    ['O3', 'O-2', '', ... ]
2781    >>> # We can also use wrapper accessors
2782    >>> atom.coordinates
2783    (0.33, 0.15, 0.5)
2784    >>> atom.refinement_flags
2785    u'FX'
2786    >>> atom.ranId
2787    4615973324315876477
2788    >>> atom.occupancy
2789    1.0
2790    """
2791    def __init__(self, data, indices, proj):
2792        self.data = data
2793        self.cx, self.ct, self.cs, self.cia = indices
2794        self.proj = proj
2795
2796    @property
2797    def label(self):
2798        '''Get the associated atom's label
2799        '''
2800        return self.data[self.ct-1]
2801
2802    @property
2803    def type(self):
2804        '''Get the associated atom's type
2805        '''
2806        return self.data[self.ct]
2807
2808    @property
2809    def refinement_flags(self):
2810        '''Get or set refinement flags for the associated atom
2811        '''
2812        return self.data[self.ct+1]
2813
2814    @refinement_flags.setter
2815    def refinement_flags(self, other):
2816        # Automatically check it is a valid refinement
2817        for c in other:
2818            if c not in ' FXU':
2819                raise ValueError("Invalid atom refinement: ", other)
2820        self.data[self.ct+1] = other
2821
2822    @property
2823    def coordinates(self):
2824        '''Get the associated atom's coordinates
2825        '''
2826        return tuple(self.data[self.cx:self.cx+3])
2827
2828    @property
2829    def occupancy(self):
2830        '''Get or set the associated atom's occupancy fraction
2831        '''
2832        return self.data[self.cx+3]
2833
2834    @occupancy.setter
2835    def occupancy(self, val):
2836        self.data[self.cx+3] = float(val)
2837
2838    @property
2839    def ranId(self):
2840        '''Get the associated atom's Random Id number
2841        '''
2842        return self.data[self.cia+8]
2843
2844    @property
2845    def adp_flag(self):
2846        '''Get the associated atom's iso/aniso setting, 'I' or 'A'
2847        '''
2848        # Either 'I' or 'A'
2849        return self.data[self.cia]
2850
2851    @property
2852    def uiso(self):
2853        '''Get or set the associated atom's Uiso or Uaniso value(s)
2854        '''
2855        if self.adp_flag == 'I':
2856            return self.data[self.cia+1]
2857        else:
2858            return self.data[self.cia+2:self.cia+8]
2859
2860    @uiso.setter
2861    def uiso(self, value):
2862        if self.adp_flag == 'I':
2863            self.data[self.cia+1] = float(value)
2864        else:
2865            assert len(value) == 6
2866            self.data[self.cia+2:self.cia+8] = [float(v) for v in value]
2867
2868class G2PwdrData(G2ObjectWrapper):
2869    """Wraps a Powder Data Histogram."""
2870    def __init__(self, data, proj, name):
2871        self.data = data
2872        self.proj = proj
2873        self.name = name
2874
2875    @staticmethod
2876    def is_valid_refinement_key(key):
2877        valid_keys = ['Limits', 'Sample Parameters', 'Background',
2878                      'Instrument Parameters']
2879        return key in valid_keys
2880
2881    #@property
2882    #def name(self):
2883    #    return self.data['data'][-1]
2884
2885    @property
2886    def ranId(self):
2887        return self.data['data'][0]['ranId']
2888
2889    @property
2890    def residuals(self):
2891        '''Provides a dictionary with with the R-factors for this histogram.
2892        Includes the weighted and unweighted profile terms (R, Rb, wR, wRb, wRmin)
2893        as well as the Bragg R-values for each phase (ph:H:Rf and ph:H:Rf^2).
2894        '''
2895        data = self.data['data'][0]
2896        return {key: data[key] for key in data
2897                if key in ['R', 'Rb', 'wR', 'wRb', 'wRmin']
2898                   or ':' in key}
2899
2900    @property
2901    def InstrumentParameters(self):
2902        '''Provides a dictionary with with the Instrument Parameters
2903        for this histogram.
2904        '''
2905        return self.data['Instrument Parameters'][0]
2906
2907    @property
2908    def SampleParameters(self):
2909        '''Provides a dictionary with with the Sample Parameters
2910        for this histogram.
2911        '''
2912        return self.data['Sample Parameters']
2913
2914    @property
2915    def Background(self):
2916        '''Provides a list with with the Background parameters
2917        for this histogram.
2918
2919        :returns: list containing a list and dict with background values
2920        '''
2921        return self.data['Background']
2922
2923    def add_back_peak(self,pos,int,sig,gam,refflags=[]):
2924        '''Adds a background peak to the Background parameters
2925       
2926        :param float pos: position of peak, a 2theta or TOF value
2927        :param float int: integrated intensity of background peak, usually large
2928        :param float sig: Gaussian width of background peak, usually large
2929        :param float gam: Lorentzian width of background peak, usually unused (small)
2930        :param list refflags: a list of 1 to 4 boolean refinement flags for
2931            pos,int,sig & gam, respectively (use [0,1] to refine int only).
2932            Defaults to [] which means nothing is refined.
2933        '''
2934        if 'peaksList' not in self.Background[1]:
2935            self.Background[1]['peaksList'] = []
2936        flags = 4*[False]
2937        for i,f in enumerate(refflags):
2938            if i>3: break
2939            flags[i] = bool(f)
2940        bpk = []
2941        for i,j in zip((pos,int,sig,gam),flags):
2942            bpk += [float(i),j]
2943        self.Background[1]['peaksList'].append(bpk)
2944        self.Background[1]['nPeaks'] = len(self.Background[1]['peaksList'])
2945
2946    def del_back_peak(self,peaknum):
2947        '''Removes a background peak from the Background parameters
2948       
2949        :param int peaknum: the number of the peak (starting from 0)
2950        '''
2951        npks = self.Background[1].get('nPeaks',0)
2952        if peaknum >= npks:
2953            raise Exception('peak {} not found in histogram {}'.format(peaknum,self.name))
2954        del self.Background[1]['peaksList'][peaknum]
2955        self.Background[1]['nPeaks'] = len(self.Background[1]['peaksList'])
2956       
2957    def ref_back_peak(self,peaknum,refflags=[]):
2958        '''Sets refinement flag for a background peak
2959       
2960        :param int peaknum: the number of the peak (starting from 0)
2961        :param list refflags: a list of 1 to 4 boolean refinement flags for
2962            pos,int,sig & gam, respectively. If a flag is not specified
2963            it defaults to False (use [0,1] to refine int only).
2964            Defaults to [] which means nothing is refined.
2965        '''
2966        npks = self.Background[1].get('nPeaks',0)
2967        if peaknum >= npks:
2968            raise Exception('peak {} not found in histogram {}'.format(peaknum,self.name))
2969        flags = 4*[False]
2970        for i,f in enumerate(refflags):
2971            if i>3: break
2972            flags[i] = bool(f)
2973        for i,f in enumerate(flags):
2974            self.Background[1]['peaksList'][peaknum][2*i+1] = f
2975                   
2976    @property
2977    def id(self):
2978        self.proj.update_ids()
2979        return self.data['data'][0]['hId']
2980
2981    @id.setter
2982    def id(self, val):
2983        self.data['data'][0]['hId'] = val
2984
2985    def fit_fixed_points(self):
2986        """Attempts to apply a background fit to the fixed points currently specified."""
2987        def SetInstParms(Inst):
2988            dataType = Inst['Type'][0]
2989            insVary = []
2990            insNames = []
2991            insVals = []
2992            for parm in Inst:
2993                insNames.append(parm)
2994                insVals.append(Inst[parm][1])
2995                if parm in ['U','V','W','X','Y','Z','SH/L','I(L2)/I(L1)','alpha',
2996                    'beta-0','beta-1','beta-q','sig-0','sig-1','sig-2','sig-q',] and Inst[parm][2]:
2997                        Inst[parm][2] = False
2998            instDict = dict(zip(insNames, insVals))
2999            instDict['X'] = max(instDict['X'], 0.01)
3000            instDict['Y'] = max(instDict['Y'], 0.01)
3001            if 'SH/L' in instDict:
3002                instDict['SH/L'] = max(instDict['SH/L'], 0.002)
3003            return dataType, instDict, insVary
3004
3005        bgrnd = self.data['Background']
3006
3007        # Need our fixed points in order
3008        bgrnd[1]['FixedPoints'].sort(key=lambda pair: pair[0])
3009        X = [x for x, y in bgrnd[1]['FixedPoints']]
3010        Y = [y for x, y in bgrnd[1]['FixedPoints']]
3011
3012        limits = self.data['Limits'][1]
3013        if X[0] > limits[0]:
3014            X = [limits[0]] + X
3015            Y = [Y[0]] + Y
3016        if X[-1] < limits[1]:
3017            X += [limits[1]]
3018            Y += [Y[-1]]
3019
3020        # Some simple lookups
3021        controls = self.proj['Controls']['data']
3022        inst, inst2 = self.data['Instrument Parameters']
3023        pwddata = self.data['data'][1]
3024
3025        # Construct the data for background fitting
3026        xBeg = np.searchsorted(pwddata[0], limits[0])
3027        xFin = np.searchsorted(pwddata[0], limits[1])
3028        xdata = pwddata[0][xBeg:xFin]
3029        ydata = si.interp1d(X,Y)(ma.getdata(xdata))
3030
3031        W = [1]*len(xdata)
3032        Z = [0]*len(xdata)
3033
3034        dataType, insDict, insVary = SetInstParms(inst)
3035        bakType, bakDict, bakVary = G2pwd.SetBackgroundParms(bgrnd)
3036
3037        # Do the fit
3038        data = np.array([xdata, ydata, W, Z, Z, Z])
3039        G2pwd.DoPeakFit('LSQ', [], bgrnd, limits, inst, inst2, data,
3040                        prevVaryList=bakVary, controls=controls)
3041
3042        # Post-fit
3043        parmDict = {}
3044        bakType, bakDict, bakVary = G2pwd.SetBackgroundParms(bgrnd)
3045        parmDict.update(bakDict)
3046        parmDict.update(insDict)
3047        pwddata[3][xBeg:xFin] *= 0
3048        pwddata[5][xBeg:xFin] *= 0
3049        pwddata[4][xBeg:xFin] = G2pwd.getBackground('', parmDict, bakType, dataType, xdata)[0]
3050
3051        # TODO adjust pwddata? GSASIIpwdGUI.py:1041
3052        # TODO update background
3053        self.proj.save()
3054
3055    def getdata(self,datatype):
3056        '''Provides access to the histogram data of the selected data type
3057
3058        :param str datatype: must be one of the following values (case is ignored)
3059       
3060           * 'X': the 2theta or TOF values for the pattern
3061           * 'Yobs': the observed intensity values
3062           * 'Yweight': the weights for each data point (1/sigma**2)
3063           * 'Ycalc': the computed intensity values
3064           * 'Background': the computed background values
3065           * 'Residual': the difference between Yobs and Ycalc (obs-calc)
3066
3067        :returns: an numpy MaskedArray with data values of the requested type
3068       
3069        '''
3070        enums = ['x', 'yobs', 'yweight', 'ycalc', 'background', 'residual']
3071        if datatype.lower() not in enums:
3072            raise G2ScriptException("Invalid datatype = "+datatype+" must be one of "+str(enums))
3073        return self.data['data'][1][enums.index(datatype.lower())]
3074       
3075    def y_calc(self):
3076        return self.data['data'][1][3]
3077
3078    def Export(self,fileroot,extension):
3079        '''Write the histogram into a file. The path is specified by fileroot and
3080        extension.
3081       
3082        :param str fileroot: name of the file, optionally with a path (extension is
3083           ignored)
3084        :param str extension: includes '.', must match an extension in global
3085           exportersByExtension['powder'] or a Exception is raised.
3086        :returns: name of file that was written
3087        '''
3088        if extension not in exportersByExtension.get('powder',[]):
3089            raise G2ScriptException('No Writer for file type = "'+extension+'"')
3090        fil = os.path.abspath(os.path.splitext(fileroot)[0]+extension)
3091        obj = exportersByExtension['powder'][extension]
3092        obj.SetFromArray(hist=self.data,histname=self.name)
3093        obj.Writer(self.name,fil)
3094        return fil
3095   
3096    def plot(self, Yobs=True, Ycalc=True, Background=True, Residual=True):
3097        try:
3098            import matplotlib.pyplot as plt
3099        except ImportError:
3100            G2fil.G2Print('Warning: Unable to import matplotlib, skipping plot')
3101            return
3102        data = self.data['data'][1]
3103        if Yobs:
3104            plt.plot(data[0], data[1], label='Yobs')
3105        if Ycalc:
3106            plt.plot(data[0], data[3], label='Ycalc')
3107        if Background:
3108            plt.plot(data[0], data[4], label='Background')
3109        if Residual:
3110            plt.plot(data[0], data[5], label="Residual")
3111
3112    def get_wR(self):
3113        """returns the overall weighted profile R factor for a histogram
3114       
3115        :returns: a wR value as a percentage or None if not defined
3116        """
3117        return self['data'][0].get('wR')
3118
3119    def set_refinements(self, refs):
3120        """Sets the histogram refinement parameter 'key' to the specification 'value'.
3121
3122        :param dict refs: A dictionary of the parameters to be set. See the
3123                          :ref:`Histogram_parameters_table` table for a description of
3124                          what these dictionaries should be.
3125
3126        :returns: None
3127
3128        """
3129        do_fit_fixed_points = False
3130        for key, value in refs.items():
3131            if key == 'Limits':
3132                old_limits = self.data['Limits'][1]
3133                new_limits = value
3134                if isinstance(new_limits, dict):
3135                    if 'low' in new_limits:
3136                        old_limits[0] = new_limits['low']
3137                    if 'high' in new_limits:
3138                        old_limits[1] = new_limits['high']
3139                else:
3140                    old_limits[0], old_limits[1] = new_limits
3141            elif key == 'Sample Parameters':
3142                sample = self.data['Sample Parameters']
3143                for sparam in value:
3144                    if sparam not in sample:
3145                        raise ValueError("Unknown refinement parameter, "
3146                                         + str(sparam))
3147                    sample[sparam][1] = True
3148            elif key == 'Background':
3149                bkg, peaks = self.data['Background']
3150
3151                # If True or False, just set the refine parameter
3152                if value in (True, False):
3153                    bkg[1] = value
3154                    return
3155
3156                if 'type' in value:
3157                    bkg[0] = value['type']
3158                if 'refine' in value:
3159                    bkg[1] = value['refine']
3160                if 'no. coeffs' in value:
3161                    cur_coeffs = bkg[2]
3162                    n_coeffs = value['no. coeffs']
3163                    if n_coeffs > cur_coeffs:
3164                        for x in range(n_coeffs - cur_coeffs):
3165                            bkg.append(0.0)
3166                    else:
3167                        for _ in range(cur_coeffs - n_coeffs):
3168                            bkg.pop()
3169                    bkg[2] = n_coeffs
3170                if 'coeffs' in value:
3171                    bkg[3:] = value['coeffs']
3172                if 'FixedPoints' in value:
3173                    peaks['FixedPoints'] = [(float(a), float(b))
3174                                            for a, b in value['FixedPoints']]
3175                if value.get('fit fixed points', False):
3176                    do_fit_fixed_points = True
3177                if 'peaks' in value:
3178                    for i,flags in enumerate(value['peaks']):
3179                        self.ref_back_peak(i,flags)
3180
3181            elif key == 'Instrument Parameters':
3182                instrument, secondary = self.data['Instrument Parameters']
3183                for iparam in value:
3184                    try:
3185                        instrument[iparam][2] = True
3186                    except IndexError:
3187                        raise ValueError("Invalid key:", iparam)
3188            else:
3189                raise ValueError("Unknown key:", key)
3190        # Fit fixed points after the fact - ensure they are after fixed points
3191        # are added
3192        if do_fit_fixed_points:
3193            # Background won't be fit if refinement flag not set
3194            orig = self.data['Background'][0][1]
3195            self.data['Background'][0][1] = True
3196            self.fit_fixed_points()
3197            # Restore the previous value
3198            self.data['Background'][0][1] = orig
3199
3200    def clear_refinements(self, refs):
3201        """Clears the refinement parameter 'key' and its associated value.
3202
3203        :param dict refs: A dictionary of parameters to clear.
3204          See the :ref:`Histogram_parameters_table` table for what can be specified.
3205        """
3206        for key, value in refs.items():
3207            if key == 'Limits':
3208                old_limits, cur_limits = self.data['Limits']
3209                cur_limits[0], cur_limits[1] = old_limits
3210            elif key == 'Sample Parameters':
3211                sample = self.data['Sample Parameters']
3212                for sparam in value:
3213                    sample[sparam][1] = False
3214            elif key == 'Background':
3215                bkg, peaks = self.data['Background']
3216
3217                # If True or False, just set the refine parameter
3218                if value in (True, False):
3219                    bkg[1] = False
3220                    return
3221
3222                bkg[1] = False
3223                if 'FixedPoints' in value:
3224                    if 'FixedPoints' in peaks:
3225                        del peaks['FixedPoints']
3226                if 'peaks' in value:
3227                    for i in range(len(self.Background[1]['peaksList'])):
3228                        self.ref_back_peak(i,[])
3229            elif key == 'Instrument Parameters':
3230                instrument, secondary = self.data['Instrument Parameters']
3231                for iparam in value:
3232                    instrument[iparam][2] = False
3233            else:
3234                raise ValueError("Unknown key:", key)
3235
3236    def add_peak(self,area,dspace=None,Q=None,ttheta=None):
3237        '''Adds a single peak to the peak list
3238        :param float area: peak area
3239        :param float dspace: peak position as d-space (A)
3240        :param float Q: peak position as Q (A-1)
3241        :param float ttheta: peak position as 2Theta (deg)
3242
3243        Note: only one of the parameters dspace, Q or ttheta may be specified
3244        '''
3245        import GSASIIlattice as G2lat
3246        import GSASIImath as G2mth
3247        if (not dspace) + (not Q) + (not ttheta) != 2:
3248            G2fil.G2Print('add_peak error: too many or no peak position(s) specified')
3249            return
3250        pos = ttheta
3251        Parms,Parms2 = self.data['Instrument Parameters']
3252        if Q:
3253            pos = G2lat.Dsp2pos(Parms,2.*np.pi/Q)
3254        elif dspace:
3255            pos = G2lat.Dsp2pos(Parms,dspace)
3256        peaks = self.data['Peak List']
3257        peaks['sigDict'] = {}        #no longer valid
3258        peaks['peaks'].append(G2mth.setPeakparms(Parms,Parms2,pos,area))
3259
3260    def set_peakFlags(self,peaklist=None,area=None,pos=None,sig=None,gam=None):
3261        '''Set refinement flags for peaks
3262       
3263        :param list peaklist: a list of peaks to change flags. If None (default), changes
3264          are made to all peaks.
3265        :param bool area: Sets or clears the refinement flag for the peak area value.
3266          If None (the default), no change is made.
3267        :param bool pos: Sets or clears the refinement flag for the peak position value.
3268          If None (the default), no change is made.
3269        :param bool sig: Sets or clears the refinement flag for the peak sig (Gaussian width) value.
3270          If None (the default), no change is made.
3271        :param bool gam: Sets or clears the refinement flag for the peak sig (Lorentzian width) value.
3272          If None (the default), no change is made.
3273         
3274        Note that when peaks are first created the area flag is on and the other flags are
3275        initially off.
3276
3277        Example::
3278       
3279           set_peakFlags(sig=False,gam=True)
3280
3281        causes the sig refinement flag to be cleared and the gam flag to be set, in both cases for
3282        all peaks. The position and area flags are not changed from their previous values.
3283        '''
3284        peaks = self.data['Peak List']
3285        if peaklist is None:
3286            peaklist = range(len(peaks['peaks']))
3287        for i in peaklist:
3288            for var,j in [(area,3),(pos,1),(sig,5),(gam,7)]:
3289                if var is not None:
3290                    peaks['peaks'][i][j] = var
3291           
3292    def refine_peaks(self):
3293        '''Causes a refinement of peak position, background and instrument parameters
3294        '''
3295        import GSASIIpwd as G2pwd
3296        controls = self.proj.data.get('Controls',{})
3297        controls = controls.get('data',
3298                            {'deriv type':'analytic','min dM/M':0.001,}     #fill in defaults if needed
3299                            )
3300        peaks = self.data['Peak List']
3301        Parms,Parms2 = self.data['Instrument Parameters']
3302        background = self.data['Background']
3303        limits = self.data['Limits'][1]
3304        bxye = np.zeros(len(self.data['data'][1][1]))
3305        peaks['sigDict'] = G2pwd.DoPeakFit('LSQ',peaks['peaks'],background,limits,
3306                                           Parms,Parms2,self.data['data'][1],bxye,[],
3307                                           False,controls,None)[0]
3308
3309    def SaveProfile(self,filename):
3310        '''Writes a GSAS-II (new style) .instprm file
3311        '''
3312        data,Parms2 = self.data['Instrument Parameters']
3313        filename = os.path.splitext(filename)[0]+'.instprm'         # make sure extension is .instprm
3314        File = open(filename,'w')
3315        File.write("#GSAS-II instrument parameter file; do not add/delete items!\n")
3316        for item in data:
3317            File.write(item+':'+str(data[item][1])+'\n')
3318        File.close()
3319        G2fil.G2Print ('Instrument parameters saved to: '+filename)
3320
3321    def LoadProfile(self,filename,bank=0):
3322        '''Reads a GSAS-II (new style) .instprm file and overwrites the current
3323        parameters
3324
3325        :param str filename: instrument parameter file name, extension ignored if not
3326          .instprm
3327        :param int bank: bank number to read, defaults to zero
3328        '''
3329        filename = os.path.splitext(filename)[0]+'.instprm'         # make sure extension is .instprm
3330        File = open(filename,'r')
3331        S = File.readline()
3332        newItems = []
3333        newVals = []
3334        Found = False
3335        while S:
3336            if S[0] == '#':
3337                if Found:
3338                    break
3339                if 'Bank' in S:
3340                    if bank == int(S.split(':')[0].split()[1]):
3341                        S = File.readline()
3342                        continue
3343                    else:
3344                        S = File.readline()
3345                        while S and '#Bank' not in S:
3346                            S = File.readline()
3347                        continue
3348                else:   #a non #Bank file
3349                    S = File.readline()
3350                    continue
3351            Found = True
3352            [item,val] = S[:-1].split(':')
3353            newItems.append(item)
3354            try:
3355                newVals.append(float(val))
3356            except ValueError:
3357                newVals.append(val)                       
3358            S = File.readline()               
3359        File.close()
3360        LoadG2fil()
3361        self.data['Instrument Parameters'][0] = G2fil.makeInstDict(newItems,newVals,len(newVals)*[False,])
3362
3363       
3364class G2Phase(G2ObjectWrapper):
3365    """A wrapper object around a given phase.
3366
3367    Author: Jackson O'Donnell (jacksonhodonnell .at. gmail.com)
3368    """
3369    def __init__(self, data, name, proj):
3370        self.data = data
3371        self.name = name
3372        self.proj = proj
3373
3374    @staticmethod
3375    def is_valid_refinement_key(key):
3376        valid_keys = ["Cell", "Atoms", "LeBail"]
3377        return key in valid_keys
3378
3379    @staticmethod
3380    def is_valid_HAP_refinement_key(key):
3381        valid_keys = ["Babinet", "Extinction", "HStrain", "Mustrain",
3382                      "Pref.Ori.", "Show", "Size", "Use", "Scale"]
3383        return key in valid_keys
3384
3385    def atom(self, atomlabel):
3386        """Returns the atom specified by atomlabel, or None if it does not
3387        exist.
3388
3389        :param str atomlabel: The name of the atom (e.g. "O2")
3390        :returns: A :class:`G2AtomRecord` object
3391            representing the atom.
3392        """
3393        # Consult GSASIIobj.py for the meaning of this
3394        cx, ct, cs, cia = self.data['General']['AtomPtrs']
3395        ptrs = [cx, ct, cs, cia]
3396        atoms = self.data['Atoms']
3397        for atom in atoms:
3398            if atom[ct-1] == atomlabel:
3399                return G2AtomRecord(atom, ptrs, self.proj)
3400
3401    def atoms(self):
3402        """Returns a list of atoms present in the current phase.
3403
3404        :returns: A list of :class:`G2AtomRecord` objects.
3405
3406        .. seealso::
3407            :meth:`~G2Phase.atom`
3408            :class:`G2AtomRecord`
3409        """
3410        ptrs = self.data['General']['AtomPtrs']
3411        output = []
3412        atoms = self.data['Atoms']
3413        for atom in atoms:
3414            output.append(G2AtomRecord(atom, ptrs, self.proj))
3415        return output
3416
3417    def histograms(self):
3418        '''Returns a list of histogram names associated with the current phase ordered
3419        as they appear in the tree (see :meth:`G2Project.histograms`).
3420        '''
3421        return [i.name for i in self.proj.histograms() if i.name in self.data.get('Histograms', {})] 
3422
3423    @property
3424    def ranId(self):
3425        return self.data['ranId']
3426
3427    @property
3428    def id(self):
3429        return self.data['pId']
3430
3431    @id.setter
3432    def id(self, val):
3433        self.data['pId'] = val
3434
3435    def get_cell(self):
3436        """Returns a dictionary of the cell parameters, with keys:
3437            'length_a', 'length_b', 'length_c', 'angle_alpha', 'angle_beta', 'angle_gamma', 'volume'
3438
3439        :returns: a dict
3440
3441        .. seealso::
3442           :meth:`~G2Phase.get_cell_and_esd`
3443
3444        """
3445        cell = self.data['General']['Cell']
3446        return {'length_a': cell[1], 'length_b': cell[2], 'length_c': cell[3],
3447                'angle_alpha': cell[4], 'angle_beta': cell[5], 'angle_gamma': cell[6],
3448                'volume': cell[7]}
3449
3450    def get_cell_and_esd(self):
3451        """
3452        Returns a pair of dictionaries, the first representing the unit cell, the second
3453        representing the estimated standard deviations of the unit cell.
3454
3455        :returns: a tuple of two dictionaries
3456
3457        .. seealso::
3458           :meth:`~G2Phase.get_cell`
3459
3460        """
3461        # translated from GSASIIstrIO.ExportBaseclass.GetCell
3462        import GSASIIlattice as G2lat
3463        import GSASIImapvars as G2mv
3464        try:
3465            pfx = str(self.id) + '::'
3466            sgdata = self['General']['SGData']
3467            covDict = self.proj['Covariance']['data']
3468
3469            parmDict = dict(zip(covDict.get('varyList',[]),
3470                                covDict.get('variables',[])))
3471            sigDict = dict(zip(covDict.get('varyList',[]),
3472                               covDict.get('sig',[])))
3473
3474            if covDict.get('covMatrix') is not None:
3475                sigDict.update(G2mv.ComputeDepESD(covDict['covMatrix'],
3476                                                  covDict['varyList'],
3477                                                  parmDict))
3478
3479            A, sigA = G2strIO.cellFill(pfx, sgdata, parmDict, sigDict)
3480            cellSig = G2strIO.getCellEsd(pfx, sgdata, A, self.proj['Covariance']['data'])
3481            cellList = G2lat.A2cell(A) + (G2lat.calc_V(A),)
3482            cellDict, cellSigDict = {}, {}
3483            for i, key in enumerate(['length_a', 'length_b', 'length_c',
3484                                     'angle_alpha', 'angle_beta', 'angle_gamma',
3485                                     'volume']):
3486                cellDict[key] = cellList[i]
3487                cellSigDict[key] = cellSig[i]
3488            return cellDict, cellSigDict
3489        except KeyError:
3490            cell = self.get_cell()
3491            return cell, {key: 0.0 for key in cell}
3492
3493    def export_CIF(self, outputname, quickmode=True):
3494        """Write this phase to a .cif file named outputname
3495
3496        :param str outputname: The name of the .cif file to write to
3497        :param bool quickmode: Currently ignored. Carryover from exports.G2export_CIF"""
3498        # This code is all taken from exports/G2export_CIF.py
3499        # Functions copied have the same names
3500        import GSASIImath as G2mth
3501        import GSASIImapvars as G2mv
3502        from exports import G2export_CIF as cif
3503
3504        CIFdate = dt.datetime.strftime(dt.datetime.now(),"%Y-%m-%dT%H:%M")
3505        CIFname = os.path.splitext(self.proj.filename)[0]
3506        CIFname = os.path.split(CIFname)[1]
3507        CIFname = ''.join([c if ord(c) < 128 else ''
3508                           for c in CIFname.replace(' ', '_')])
3509        try:
3510            author = self.proj['Controls']['data'].get('Author','').strip()
3511        except KeyError:
3512            pass
3513        oneblock = True
3514
3515        covDict = self.proj['Covariance']['data']
3516        parmDict = dict(zip(covDict.get('varyList',[]),
3517                            covDict.get('variables',[])))
3518        sigDict = dict(zip(covDict.get('varyList',[]),
3519                           covDict.get('sig',[])))
3520
3521        if covDict.get('covMatrix') is not None:
3522            sigDict.update(G2mv.ComputeDepESD(covDict['covMatrix'],
3523                                              covDict['varyList'],
3524                                              parmDict))
3525
3526        with open(outputname, 'w') as fp:
3527            fp.write(' \n' + 70*'#' + '\n')
3528            cif.WriteCIFitem(fp, 'data_' + CIFname)
3529            # from exports.G2export_CIF.WritePhaseInfo
3530            cif.WriteCIFitem(fp, '\n# phase info for '+str(self.name) + ' follows')
3531            cif.WriteCIFitem(fp, '_pd_phase_name', self.name)
3532            # TODO get esds
3533            cellDict = self.get_cell()
3534            defsigL = 3*[-0.00001] + 3*[-0.001] + [-0.01] # significance to use when no sigma
3535            names = ['length_a','length_b','length_c',
3536                     'angle_alpha','angle_beta ','angle_gamma',
3537                     'volume']
3538            for key, val in cellDict.items():
3539                cif.WriteCIFitem(fp, '_cell_' + key, G2mth.ValEsd(val))
3540
3541            cif.WriteCIFitem(fp, '_symmetry_cell_setting',
3542                         self.data['General']['SGData']['SGSys'])
3543
3544            spacegroup = self.data['General']['SGData']['SpGrp'].strip()
3545            # regularize capitalization and remove trailing H/R
3546            spacegroup = spacegroup[0].upper() + spacegroup[1:].lower().rstrip('rh ')
3547            cif.WriteCIFitem(fp, '_symmetry_space_group_name_H-M', spacegroup)
3548
3549            # generate symmetry operations including centering and center of symmetry
3550            SymOpList, offsetList, symOpList, G2oprList, G2opcodes = G2spc.AllOps(
3551                self.data['General']['SGData'])
3552            cif.WriteCIFitem(fp, 'loop_\n    _space_group_symop_id\n    _space_group_symop_operation_xyz')
3553            for i, op in enumerate(SymOpList,start=1):
3554                cif.WriteCIFitem(fp, '   {:3d}  {:}'.format(i,op.lower()))
3555
3556            # TODO skipped histograms, exports/G2export_CIF.py:880
3557
3558            # report atom params
3559            if self.data['General']['Type'] in ['nuclear','macromolecular']:        #this needs macromolecular variant, etc!
3560                cif.WriteAtomsNuclear(fp, self.data, self.name, parmDict, sigDict, [])
3561                # self._WriteAtomsNuclear(fp, parmDict, sigDict)
3562            else:
3563                raise G2ScriptException("no export for "+str(self.data['General']['Type'])+" coordinates implemented")
3564            # report cell contents
3565            cif.WriteComposition(fp, self.data, self.name, parmDict)
3566            if not quickmode and self.data['General']['Type'] == 'nuclear':      # report distances and angles
3567                # WriteDistances(fp,self.name,SymOpList,offsetList,symOpList,G2oprList)
3568                raise NotImplementedError("only quickmode currently supported")
3569            if 'Map' in self.data['General'] and 'minmax' in self.data['General']['Map']:
3570                cif.WriteCIFitem(fp,'\n# Difference density results')
3571                MinMax = self.data['General']['Map']['minmax']
3572                cif.WriteCIFitem(fp,'_refine_diff_density_max',G2mth.ValEsd(MinMax[0],-0.009))
3573                cif.WriteCIFitem(fp,'_refine_diff_density_min',G2mth.ValEsd(MinMax[1],-0.009))
3574
3575
3576    def set_refinements(self, refs):
3577        """Sets the phase refinement parameter 'key' to the specification 'value'
3578
3579        :param dict refs: A dictionary of the parameters to be set. See the
3580                          :ref:`Phase_parameters_table` table for a description of
3581                          this dictionary.
3582
3583        :returns: None"""
3584        for key, value in refs.items():
3585            if key == "Cell":
3586                self.data['General']['Cell'][0] = value
3587
3588            elif key == "Atoms":
3589                for atomlabel, atomrefinement in value.items():
3590                    if atomlabel == 'all':
3591                        for atom in self.atoms():
3592                            atom.refinement_flags = atomrefinement
3593                    else:
3594                        atom = self.atom(atomlabel)
3595                        if atom is None:
3596                            raise ValueError("No such atom: " + atomlabel)
3597                        atom.refinement_flags = atomrefinement
3598
3599            elif key == "LeBail":
3600                hists = self.data['Histograms']
3601                for hname, hoptions in hists.items():
3602                    if 'LeBail' not in hoptions:
3603                        hoptions['newLeBail'] = bool(True)
3604                    hoptions['LeBail'] = bool(value)
3605            else:
3606                raise ValueError("Unknown key:", key)
3607
3608    def clear_refinements(self, refs):
3609        """Clears a given set of parameters.
3610
3611        :param dict refs: The parameters to clear.
3612          See the :ref:`Phase_parameters_table` table for what can be specified.
3613        """
3614        for key, value in refs.items():
3615            if key == "Cell":
3616                self.data['General']['Cell'][0] = False
3617            elif key == "Atoms":
3618                cx, ct, cs, cia = self.data['General']['AtomPtrs']
3619
3620                for atomlabel in value:
3621                    atom = self.atom(atomlabel)
3622                    # Set refinement to none
3623                    atom.refinement_flags = ' '
3624            elif key == "LeBail":
3625                hists = self.data['Histograms']
3626                for hname, hoptions in hists.items():
3627                    if 'LeBail' not in hoptions:
3628                        hoptions['newLeBail'] = True
3629                    hoptions['LeBail'] = False
3630            else:
3631                raise ValueError("Unknown key:", key)
3632
3633    def set_HAP_refinements(self, refs, histograms='all'):
3634        """Sets the given HAP refinement parameters between the current phase and
3635        the specified histograms.
3636
3637        :param dict refs: A dictionary of the parameters to be set. See
3638                          the :ref:`HAP_parameters_table` table for a description of this
3639                          dictionary.
3640        :param histograms: Either 'all' (default) or a list of the histograms by index, name
3641            or object. The index number is relative to all histograms in the tree, not to
3642            those in the phase.
3643            Histograms not associated with the current phase will be ignored.
3644            whose HAP parameters will be set with this phase. Histogram and phase
3645            must already be associated.
3646        :returns: None
3647        """
3648        if not self.data.get('Histograms',[]):
3649            G2fil.G2Print("Error likely: Phase {} has no linked histograms".format(self.name))
3650            return
3651           
3652        if histograms == 'all':
3653            histograms = self.data['Histograms'].keys()
3654        else:
3655            histograms = [self._decodeHist(h) for h in histograms
3656                          if self._decodeHist(h) in self.data['Histograms']]   
3657        if not histograms:
3658            G2fil.G2Print("Warning: Skipping HAP set for phase {}, no selected histograms".format(self.name))
3659            return
3660        for key, val in refs.items():
3661            if key == 'Babinet':
3662                try:
3663                    sets = list(val)
3664                except ValueError:
3665                    sets = ['BabA', 'BabU']
3666                for param in sets:
3667                    if param not in ['BabA', 'BabU']:
3668                        raise ValueError("Not sure what to do with" + param)
3669                    for h in histograms:
3670                        self.data['Histograms'][h]['Babinet'][param][1] = True
3671            elif key == 'Extinction':
3672                for h in histograms:
3673                    self.data['Histograms'][h]['Extinction'][1] = bool(val)
3674            elif key == 'HStrain':
3675                if isinstance(val,list) or isinstance(val,tuple):
3676                    for h in histograms:
3677                        if len(self.data['Histograms'][h]['HStrain'][1]) != len(val):
3678                            raise Exception('Need {} HStrain terms for phase {} hist {}'
3679                                .format(len(self.data['Histograms'][h]['HStrain'][1]),self.name,h))
3680                        for i,v in enumerate(val):
3681                            self.data['Histograms'][h]['HStrain'][1][i] = bool(v)
3682                else:
3683                    for h in histograms:
3684                        self.data['Histograms'][h]['HStrain'][1] = [bool(val) for p in self.data['Histograms'][h]['HStrain'][1]]
3685            elif key == 'Mustrain':
3686                for h in histograms:
3687                    mustrain = self.data['Histograms'][h]['Mustrain']
3688                    newType = None
3689                    direction = None
3690                    if isinstance(val, strtypes):
3691                        if val in ['isotropic', 'uniaxial', 'generalized']:
3692                            newType = val
3693                        else:
3694                            raise ValueError("Not a Mustrain type: " + val)
3695                    elif isinstance(val, dict):
3696                        newType = val.get('type', None)
3697                        direction = val.get('direction', None)
3698
3699                    if newType:
3700                        mustrain[0] = newType
3701                        if newType == 'isotropic':
3702                            mustrain[2][0] = True == val.get('refine',False)
3703                            mustrain[5] = [False for p in mustrain[4]]
3704                        elif newType == 'uniaxial':
3705                            if 'refine' in val:
3706                                mustrain[2][0] = False
3707                                types = val['refine']
3708                                if isinstance(types, strtypes):
3709                                    types = [types]
3710                                elif isinstance(types, bool):
3711                                    mustrain[2][1] = types
3712                                    mustrain[2][2] = types
3713                                    types = []
3714                                else:
3715                                    raise ValueError("Not sure what to do with: "
3716                                                     + str(types))
3717                            else:
3718                                types = []
3719
3720                            for unitype in types:
3721                                if unitype == 'equatorial':
3722                                    mustrain[2][0] = True
3723                                elif unitype == 'axial':
3724                                    mustrain[2][1] = True
3725                                else:
3726                                    msg = 'Invalid uniaxial mustrain type'
3727                                    raise ValueError(msg + ': ' + unitype)
3728                        else:  # newtype == 'generalized'
3729                            mustrain[2] = [False for p in mustrain[1]]
3730                            if 'refine' in val:
3731                                mustrain[5] = [True == val['refine']]*len(mustrain[5])
3732
3733                    if direction:
3734                        if len(direction) != 3:
3735                            raise ValueError("Expected hkl, found", direction)
3736                        direction = [int(n) for n in direction]
3737                        mustrain[3] = direction
3738            elif key == 'Size':
3739                newSize = None
3740                if 'value' in val:
3741                    newSize = float(val['value'])
3742                for h in histograms:
3743                    size = self.data['Histograms'][h]['Size']
3744                    newType = None
3745                    direction = None
3746                    if isinstance(val, strtypes):
3747                        if val in ['isotropic', 'uniaxial', 'ellipsoidal']:
3748                            newType = val
3749                        else:
3750                            raise ValueError("Not a valid Size type: " + val)
3751                    elif isinstance(val, dict):
3752                        newType = val.get('type', None)
3753                        direction = val.get('direction', None)
3754
3755                    if newType:
3756                        size[0] = newType
3757                        refine = bool(val.get('refine'))
3758                        if newType == 'isotropic' and refine is not None:
3759                            size[2][0] = bool(refine)
3760                            if newSize: size[1][0] = newSize
3761                        elif newType == 'uniaxial' and refine is not None:
3762                            size[2][1] = bool(refine)
3763                            size[2][2] = bool(refine)
3764                            if newSize: size[1][1] = size[1][2] =newSize
3765                        elif newType == 'ellipsoidal' and refine is not None:
3766                            size[5] = [bool(refine) for p in size[5]]
3767                            if newSize: size[4] = [newSize for p in size[4]]
3768
3769                    if direction:
3770                        if len(direction) != 3:
3771                            raise ValueError("Expected hkl, found", direction)
3772                        direction = [int(n) for n in direction]
3773                        size[3] = direction
3774            elif key == 'Pref.Ori.':
3775                for h in histograms:
3776                    self.data['Histograms'][h]['Pref.Ori.'][2] = bool(val)
3777            elif key == 'Show':
3778                for h in histograms:
3779                    self.data['Histograms'][h]['Show'] = bool(val)
3780            elif key == 'Use':
3781                for h in histograms:
3782                    self.data['Histograms'][h]['Use'] = bool(val)
3783            elif key == 'Scale':
3784                for h in histograms:
3785                    self.data['Histograms'][h]['Scale'][1] = bool(val)
3786            else:
3787                G2fil.G2Print(u'Warning: Unknown HAP key: '+key)
3788
3789    def clear_HAP_refinements(self, refs, histograms='all'):
3790        """Clears the given HAP refinement parameters between this phase and
3791        the given histograms.
3792
3793        :param dict refs: A dictionary of the parameters to be cleared.
3794            See the the :ref:`HAP_parameters_table` table for what can be specified.
3795        :param histograms: Either 'all' (default) or a list of the histograms by index, name
3796            or object.
3797            The index number is relative to all histograms in the tree, not to
3798            those in the phase.
3799            Histograms not associated with the current phase will be ignored.
3800            whose HAP parameters will be set with this phase. Histogram and phase
3801            must already be associated
3802        :returns: None
3803        """
3804        if histograms == 'all':
3805            histograms = self.data['Histograms'].keys()
3806        else:
3807            histograms = [self._decodeHist(h) for h in histograms
3808                          if self._decodeHist(h) in self.data['Histograms']]   
3809
3810        for key, val in refs.items():
3811            for h in histograms:
3812                if key == 'Babinet':
3813                    try:
3814                        sets = list(val)
3815                    except ValueError:
3816                        sets = ['BabA', 'BabU']
3817                    for param in sets:
3818                        if param not in ['BabA', 'BabU']:
3819                            raise ValueError("Not sure what to do with" + param)
3820                        for h in histograms:
3821                            self.data['Histograms'][h]['Babinet'][param][1] = False
3822                elif key == 'Extinction':
3823                    for h in histograms:
3824                        self.data['Histograms'][h]['Extinction'][1] = False
3825                elif key == 'HStrain':
3826                    for h in histograms:
3827                        self.data['Histograms'][h]['HStrain'][1] = [False for p in self.data['Histograms'][h]['HStrain'][1]]
3828                elif key == 'Mustrain':
3829                    for h in histograms:
3830                        mustrain = self.data['Histograms'][h]['Mustrain']
3831                        mustrain[2] = [False for p in mustrain[2]]
3832                        mustrain[5] = [False for p in mustrain[4]]
3833                elif key == 'Pref.Ori.':
3834                    for h in histograms:
3835                        self.data['Histograms'][h]['Pref.Ori.'][2] = False
3836                elif key == 'Show':
3837                    for h in histograms:
3838                        self.data['Histograms'][h]['Show'] = False
3839                elif key == 'Size':
3840                    for h in histograms:
3841                        size = self.data['Histograms'][h]['Size']
3842                        size[2] = [False for p in size[2]]
3843                        size[5] = [False for p in size[5]]
3844                elif key == 'Use':
3845                    for h in histograms:
3846                        self.data['Histograms'][h]['Use'] = False
3847                elif key == 'Scale':
3848                    for h in histograms:
3849                        self.data['Histograms'][h]['Scale'][1] = False
3850                else:
3851                    G2fil.G2Print(u'Warning: Unknown HAP key: '+key)
3852
3853    def _decodeHist(self,hist):
3854        '''Convert a histogram reference to a histogram name string
3855        '''
3856        if isinstance(hist, G2PwdrData):
3857            return hist.name
3858        elif hist in self.data['Histograms']:
3859            return hist
3860        elif type(hist) is int:
3861            return self.proj.histograms()[hist].name
3862        else:
3863            raise G2ScriptException("Invalid histogram reference: "+str(hist))
3864       
3865    def getHAPvalues(self, histname):
3866        """Returns a dict with HAP values for the selected histogram
3867
3868        :param histogram: is a histogram object (:class:`G2PwdrData`) or
3869            a histogram name or the index number of the histogram.
3870            The index number is relative to all histograms in the tree, not to
3871            those in the phase.
3872        :returns: HAP value dict
3873        """
3874        return self.data['Histograms'][self._decodeHist(histname)]
3875
3876    def copyHAPvalues(self, sourcehist, targethistlist='all', skip=[], use=None):
3877        """Copies HAP parameters for one histogram to a list of other histograms.
3878        Use skip or use to select specific entries to be copied or not used.
3879
3880        :param sourcehist: is a histogram object (:class:`G2PwdrData`) or
3881            a histogram name or the index number of the histogram to copy
3882            parameters from.
3883            The index number is relative to all histograms in the tree, not to
3884            those in the phase.
3885        :param list targethistlist: a list of histograms where each item in the
3886            list can be a histogram object (:class:`G2PwdrData`),
3887            a histogram name or the index number of the histogram.
3888            If the string 'all' (default), then all histograms in the phase
3889            are used.
3890        :param list skip: items in the HAP dict that should not be
3891            copied. The default is an empty list, which causes all items
3892            to be copied. To see a list of items in the dict, use
3893            :meth:`getHAPvalues` or use an invalid item, such as '?'.
3894        :param list use: specifies the items in the HAP dict should be
3895            copied. The default is None, which causes all items
3896            to be copied.
3897
3898        examples::
3899
3900            ph0.copyHAPvalues(0,[1,2,3])
3901            ph0.copyHAPvalues(0,use=['HStrain','Size'])
3902
3903        The first example copies all HAP parameters from the first histogram to
3904        the second, third and fourth histograms (as listed in the project tree).
3905        The second example copies only the 'HStrain' (Dij parameters and
3906        refinement flags) and the 'Size' (crystallite size settings, parameters
3907        and refinement flags) from the first histogram to all histograms.
3908        """
3909        sourcehist = self._decodeHist(sourcehist)
3910        if targethistlist == 'all':
3911            targethistlist = self.histograms()
3912       
3913        copydict = copy.deepcopy(self.data['Histograms'][sourcehist])
3914        for item in skip:
3915            if item in list(copydict.keys()):
3916                del copydict[item]
3917            else:
3918                G2fil.G2Print('items in HAP dict are: {}'.format(
3919                    list(self.data['Histograms'][sourcehist])))
3920                raise Exception('HAP skip list entry {} invalid'.format(item))
3921        if use:
3922            for item in list(copydict.keys()):
3923                if item not in use:
3924                    del copydict[item]
3925
3926        G2fil.G2Print('Copying item(s) {} from histogram {}'.format(list(copydict.keys()),sourcehist))
3927        G2fil.G2Print(' to histogram(s) {}'.format([self._decodeHist(h) for h in targethistlist]))
3928        for h in targethistlist:
3929            h = self._decodeHist(h)
3930            if h not in self.data['Histograms']:
3931                G2fil.G2Print('Unexpected Warning: histogram {} not in phase {}'.format(h,self.name))
3932                continue
3933            self.data['Histograms'][h].update(copy.deepcopy(copydict))
3934
3935    def setHAPvalues(self, HAPdict, targethistlist='all', skip=[], use=None):
3936        """Copies HAP parameters for one histogram to a list of other histograms.
3937        Use skip or use to select specific entries to be copied or not used.
3938        Note that ``HStrain`` and sometimes ``Mustrain`` values can be specific to
3939        a Laue class and should be copied with care between phases of different
3940        symmetry. A "sanity check" on the number of Dij terms is made if ``HStrain``
3941        values are copied.
3942
3943        :param dict HAPdict: is a dict returned by :meth:`getHAPvalues` containing
3944            HAP parameters.
3945        :param list targethistlist: a list of histograms where each item in the
3946            list can be a histogram object (:class:`G2PwdrData`),
3947            a histogram name or the index number of the histogram.
3948            The index number is relative to all histograms in the tree, not to
3949            those in the phase.
3950            If the string 'all' (default), then all histograms in the phase
3951            are used.
3952        :param list skip: items in the HAP dict that should not be
3953            copied. The default is an empty list, which causes all items
3954            to be copied. To see a list of items in the dict, use
3955            :meth:`getHAPvalues` or use an invalid item, such as '?'.
3956        :param list use: specifies the items in the HAP dict should be
3957            copied. The default is None, which causes all items
3958            to be copied.
3959
3960        example::
3961
3962            HAPdict = ph0.getHAPvalues(0)
3963            ph1.setHAPvalues(HAPdict,use=['HStrain','Size'])
3964
3965        This copies the Dij (hydrostatic strain) HAP parameters and the
3966        crystallite size broadening terms from the first histogram in
3967        phase ``ph0`` to all histograms in phase ``ph1``.
3968        """
3969        if targethistlist == 'all':
3970            targethistlist = self.histograms()
3971        copydict = copy.deepcopy(HAPdict)
3972        for item in skip:
3973            if item in list(copydict.keys()):
3974                del copydict[item]
3975            else:
3976                G2fil.G2Print('items in HAP dict are: {}'.format(
3977                    list(self.data['Histograms'][sourcehist])))
3978                raise Exception('HAP skip list entry {} invalid'.format(item))
3979        if use:
3980            for item in list(copydict.keys()):
3981                if item not in use:
3982                    del copydict[item]
3983
3984        first = True
3985        for h in targethistlist:
3986            h = self._decodeHist(h)
3987            if h not in self.data['Histograms']:
3988                G2fil.G2Print('Warning: histogram {} not in phase {}'.format(h,self.name))
3989                continue
3990            if first:
3991                first = False
3992                if 'HStrain' in self.data['Histograms'][h] and 'HStrain' in copydict:
3993                    if len(copydict['HStrain'][0]) != len(self.data['Histograms'][h]['HStrain'][0]):
3994                        G2fil.G2Print('Error: HStrain has differing numbers of terms. Input: {}, phase {}: {}'.
3995                                  format(len(copydict['HStrain'][0]),
3996                                        self.name,len(self.data['Histograms'][h]['HStrain'][0])))
3997                        raise Exception('HStrain has differing numbers of terms.')
3998            self.data['Histograms'][h].update(copy.deepcopy(copydict))           
3999        G2fil.G2Print('Copied item(s) {} from dict'.format(list(copydict.keys())))
4000        G2fil.G2Print(' to histogram(s) {}'.format([self._decodeHist(h) for h in targethistlist]))
4001
4002class G2SeqRefRes(G2ObjectWrapper):
4003    '''Wrapper for a Sequential Refinement Results tree entry, containing the
4004    results for a refinement
4005
4006    As an example::
4007
4008        from __future__ import division, print_function
4009        import os,sys
4010        sys.path.insert(0,'/Users/toby/software/G2/GSASII')
4011        PathWrap = lambda fil: os.path.join('/Users/toby/Scratch/SeqTut2019Mar',fil)
4012        import GSASIIscriptable as G2sc
4013        gpx = G2sc.G2Project(PathWrap('scr4.gpx'))
4014        seq = gpx.seqref()
4015        lbl = ('a','b','c','alpha','beta','gamma','Volume')
4016        for j,h in enumerate(seq.histograms()):
4017            cell,cellU,uniq = seq.get_cell_and_esd(1,h)
4018            print(h)
4019            print([cell[i] for i in list(uniq)+[6]])
4020            print([cellU[i] for i in list(uniq)+[6]])
4021            print('')
4022        print('printed',[lbl[i] for i in list(uniq)+[6]])
4023
4024    .. seealso::
4025        :meth:`G2Project.seqref`
4026    '''
4027    def __init__(self, data, proj):
4028        self.data = data
4029        self.proj = proj
4030        self.newCellDict = {} # dict with recp. cell tensor & Dij name
4031        # newAtomDict = {} # dict with atom positions; relative & absolute
4032        for name in self.data['histNames']:
4033            self.newCellDict.update(self.data[name].get('newCellDict',{}))
4034            #newAtomDict.update(self.data[name].get('newAtomDict',{}) # dict with atom positions; relative & absolute
4035        #ESDlookup = {self.newCellDict[item][0]:item for item in self.newCellDict}
4036        #Dlookup = {item:self.newCellDict[item][0] for item in self.newCellDict}
4037        # Possible error: the next might need to be data[histNames[0]]['varyList']
4038        #atomLookup = {newAtomDict[item][0]:item for item in newAtomDict if item in self.data['varyList']}
4039        #Dlookup.update({atomLookup[parm]:parm for parm in atomLookup}
4040        #ESDlookup.update({parm:atomLookup[parm] for parm in atomLookup})
4041
4042#    @property
4043    def histograms(self):
4044        '''returns a list of histograms in the squential fit
4045        '''
4046        return self.data['histNames']
4047
4048    def get_cell_and_esd(self,phase,hist):
4049        '''Returns a vector of cell lengths and esd values
4050
4051        :param phase: A phase, which may be specified as a phase object
4052          (see :class:`G2Phase`), the phase name (str) or the index number (int)
4053          of the phase in the project, numbered starting from 0.
4054        :param hist: Specify a histogram or using the histogram name (str)
4055          or the index number (int) of the histogram in the sequential
4056          refinement (not the project), numbered as in in the project tree
4057          starting from 0.
4058        :returns: cell,cellESD,uniqCellIndx where cell (list)
4059          with the unit cell parameters (a,b,c,alpha,beta,gamma,Volume);
4060          cellESD are the standard uncertainties on the 7 unit cell
4061          parameters; and uniqCellIndx is a tuple with indicies for the
4062          unique (non-symmetry determined) unit parameters (e.g.
4063          [0,2] for a,c in a tetragonal cell)
4064        '''
4065        def striphist(var,insChar=''):
4066            'strip a histogram number from a var name'
4067            sv = var.split(':')
4068            if len(sv) <= 1: return var
4069            if sv[1]:
4070                sv[1] = insChar
4071            return ':'.join(sv)
4072        import GSASIIlattice as G2lat
4073        import GSASIIstrIO as G2stIO
4074       
4075        uniqCellLookup = [
4076        [['m3','m3m'],(0,)],
4077        [['3R','3mR'],(0,3)],
4078        [['3','3m1','31m','6/m','6/mmm','4/m','4/mmm'],(0,2)],
4079        [['mmm'],(0,1,2)],
4080        [['2/m'+'a'],(0,1,2,3)],
4081        [['2/m'+'b'],(0,1,2,4)],
4082        [['2/m'+'c'],(0,1,2,5)],
4083        [['-1'],(0,1,2,3,4,5)],
4084        ]
4085
4086        seqData,histData = self.RefData(hist)
4087        hId = histData['data'][0]['hId']
4088        phasedict = self.proj.phase(phase).data
4089        pId = phasedict['pId']
4090        pfx = str(pId)+'::' # prefix for A values from phase
4091        phfx = '%d:%d:'%(pId,hId) # Dij prefix
4092        # get unit cell & symmetry for phase
4093        RecpCellTerms = G2lat.cell2A(phasedict['General']['Cell'][1:7])
4094        zeroDict = {pfx+'A'+str(i):0.0 for i in range(6)}
4095        SGdata = phasedict['General']['SGData']
4096        # determine the cell items not defined by symmetry
4097        laue = SGdata['SGLaue'][:]
4098        if laue == '2/m':
4099            laue += SGdata['SGUniq']
4100        for symlist,celllist in uniqCellLookup:
4101            if laue in symlist:
4102                uniqCellIndx = celllist
4103                break
4104        else: # should not happen
4105            uniqCellIndx = list(range(6))
4106        initialCell = {}
4107        for i in uniqCellIndx:
4108            initialCell[str(pId)+'::A'+str(i)] =  RecpCellTerms[i]
4109
4110        esdLookUp = {}
4111        dLookup = {}
4112        # note that varyList keys are p:h:Dij while newCellDict keys are p::Dij
4113        for nKey in seqData['newCellDict']:
4114            p = nKey.find('::')+1
4115            vKey = nKey[:p] + str(hId) + nKey[p:]
4116            if vKey in seqData['varyList']:
4117                esdLookUp[self.newCellDict[nKey][0]] = nKey
4118                dLookup[nKey] = self.newCellDict[nKey][0]
4119        covData = {'varyList': [dLookup.get(striphist(v),v) for v in seqData['varyList']],
4120                'covMatrix': seqData['covMatrix']}
4121        A = RecpCellTerms[:] # make copy of starting A values
4122
4123        for i,j in enumerate(('D11','D22','D33','D12','D13','D23')):
4124            var = pfx+'A'+str(i)
4125            Dvar = phfx+j
4126            # apply Dij value if non-zero
4127            if Dvar in seqData['parmDict']:
4128                A[i] += seqData['parmDict'][Dvar]
4129            # override with fit result if is Dij varied
4130            try:
4131                A[i] = seqData['newCellDict'][esdLookUp[var]][1] # get refined value
4132            except KeyError:
4133                pass
4134        Albls = [pfx+'A'+str(i) for i in range(6)]
4135        cellDict = dict(zip(Albls,A))
4136
4137       
4138        A,zeros = G2stIO.cellFill(pfx,SGdata,cellDict,zeroDict)
4139        # convert to direct cell
4140        c = G2lat.A2cell(A)
4141        vol = G2lat.calc_V(A)
4142        cE = G2stIO.getCellEsd(pfx,SGdata,A,covData)
4143        return list(c)+[vol],cE,uniqCellIndx
4144   
4145    def get_VaryList(self,hist):
4146        '''Returns a list of the refined variables in the
4147        last refinement cycle for the selected histogram
4148
4149        :param hist: Specify a histogram or using the histogram name (str)
4150          or the index number (int) of the histogram in the sequential
4151          refinement (not the project), numbered starting from 0.
4152        :returns: a list of variables or None if no refinement has been
4153          performed.
4154        '''
4155        try:
4156            seqData,histData = self.RefData(hist)
4157            return seqData['varyList']
4158        except:
4159            return
4160
4161    def get_ParmList(self,hist):
4162        '''Returns a list of all the parameters defined in the
4163        last refinement cycle for the selected histogram
4164
4165        :param hist: Specify a histogram or using the histogram name (str)
4166          or the index number (int) of the histogram in the sequential
4167          refinement (not the project), numbered as in the project tree
4168          starting from 0.
4169        :returns: a list of parameters or None if no refinement has been
4170          performed.
4171        '''
4172        try:
4173            seqData,histData = self.RefData(hist)
4174            return list(seqData['parmDict'].keys())
4175        except:
4176            return
4177
4178    def get_Variable(self,hist,var):
4179        '''Returns the value and standard uncertainty (esd) for a variable
4180        parameters, as defined for the selected histogram
4181        in the last sequential refinement cycle
4182
4183        :param hist: Specify a histogram or using the histogram name (str)
4184          or the index number (int) of the histogram in the sequential
4185          refinement (not the project), numbered as in the project tree
4186          starting from 0.
4187        :param str var: a variable name of form '<p>:<h>:<name>', such as
4188          ':0:Scale'
4189        :returns: (value,esd) if the parameter is refined or
4190          (value, None) if the variable is in a constraint or is not
4191          refined or None if the parameter is not found.
4192        '''
4193        try:
4194            seqData,histData = self.RefData(hist)
4195            val = seqData['parmDict'][var]
4196        except:
4197            return
4198        try:
4199            pos = seqData['varyList'].index(var)
4200            esd = np.sqrt(seqData['covMatrix'][pos,pos])
4201            return (val,esd)
4202        except ValueError:
4203            return (val,None)
4204
4205    def get_Covariance(self,hist,varList):
4206        '''Returns the values and covariance matrix for a series of variable
4207        parameters, as defined for the selected histogram
4208        in the last sequential refinement cycle
4209
4210        :param hist: Specify a histogram or using the histogram name (str)
4211          or the index number (int) of the histogram in the sequential
4212          refinement (not the project), numbered as in the project tree
4213          starting from 0.
4214        :param tuple varList: a list of variable names of form '<p>:<h>:<name>'
4215        :returns: (valueList,CovMatrix) where valueList contains the (n) values
4216          in the same order as varList (also length n) and CovMatrix is a
4217          (n x n) matrix. If any variable name is not found in the varyList
4218          then None is returned.
4219
4220        Use this code, where sig provides standard uncertainties for
4221        parameters and where covArray provides the correlation between
4222        off-diagonal terms::
4223
4224            sig = np.sqrt(np.diag(covMatrix))
4225            xvar = np.outer(sig,np.ones_like(sig))
4226            covArray = np.divide(np.divide(covMatrix,xvar),xvar.T)
4227
4228        '''
4229        try:
4230            seqData,histData = self.RefData(hist)
4231        except:
4232            G2fil.G2Print('Warning: Histogram {} not found in the sequential fit'.format(hist))
4233            return
4234        missing = [i for i in varList if i not in seqData['varyList']]
4235        if missing:
4236            G2fil.G2Print('Warning: Variable(s) {} were not found in the varyList'.format(missing))
4237            return None
4238        vals = [seqData['parmDict'][i] for i in varList]
4239        import GSASIImath as G2mth
4240        cov = G2mth.getVCov(varList,seqData['varyList'],seqData['covMatrix'])
4241        return (vals,cov)
4242   
4243    def RefData(self,hist):
4244        '''Provides access to the output from a particular histogram
4245
4246        :param hist: Specify a histogram or using the histogram name (str)
4247          or the index number (int) of the histogram in the sequential
4248          refinement (not the project), numbered as in the project tree
4249          starting from 0.
4250        :returns: a list of dicts where the first element has sequential
4251          refinement results and the second element has the contents of
4252          the histogram tree items.
4253        '''
4254        try:
4255            hist = self.data['histNames'][hist]
4256        except IndexError:
4257            raise Exception('Histogram #{} is out of range from the Sequential Refinement'
4258                                .format(hist))
4259        except TypeError:
4260            pass
4261        if hist not in self.data['histNames']:
4262            raise Exception('Histogram {} is not included in the Sequential Refinement'
4263                                .format(hist))
4264        return self.data[hist],self.proj.histogram(hist).data
4265
4266class G2PDF(G2ObjectWrapper):
4267    """Wrapper for a PDF tree entry, containing the information needed to
4268    compute a PDF and the S(Q), G(r) etc. after the computation is done.
4269    Note that in a GSASIIscriptable script, instances of G2PDF will be created by
4270    calls to :meth:`G2Project.add_PDF` or :meth:`G2Project.pdf`, not via calls
4271    to :meth:`G2PDF.__init__`.
4272
4273    Example use of :class:`G2PDF`::
4274
4275       gpx.add_PDF('250umSiO2.pdfprm',0)
4276       pdf.set_formula(['Si',1],['O',2])
4277       pdf.set_background('Container',1,-0.21)
4278       for i in range(5):
4279           if pdf.optimize(): break
4280       pdf.calculate()
4281       pdf.export(gpx.filename,'S(Q), pdfGUI')
4282       gpx.save('pdfcalc.gpx')
4283
4284    .. seealso::
4285        :meth:`G2Project.pdf`
4286        :meth:`G2Project.pdfs`
4287    """
4288    def __init__(self, data, name, proj):
4289        self.data = data
4290        self.name = name
4291        self.proj = proj
4292    def set_background(self,btype,histogram,mult=-1.,refine=False):
4293        '''Sets a histogram to be used as the 'Sample Background',
4294        the 'Container' or the 'Container Background.'
4295
4296        :param str btype: Type of background to set, must contain
4297          the string 'samp' for Sample Background', 'cont' and 'back'
4298          for the 'Container Background' or only 'cont' for the
4299          'Container'. Note that capitalization and extra characters
4300          are ignored, so the full strings (such as 'Sample
4301          Background' & 'Container Background') can be used.
4302       
4303        :param histogram: A reference to a histogram,
4304          which can be reference by object, name, or number.
4305        :param float mult: a multiplier for the histogram; defaults
4306          to -1.0
4307        :param bool refine: a flag to enable refinement (only
4308          implemented for 'Sample Background'); defaults to False
4309        '''
4310        if 'samp' in btype.lower():
4311            key = 'Sample Bkg.'
4312        elif 'cont' in btype.lower() and 'back' in btype.lower():
4313            key = 'Container Bkg.'
4314        elif 'cont' in btype.lower():
4315            key = 'Container'
4316        else:
4317            raise Exception('btype = {} is invalid'.format(btype))
4318        self.data['PDF Controls'][key]['Name'] = self.proj.histogram(histogram).name
4319        self.data['PDF Controls'][key]['Mult'] = mult
4320        self.data['PDF Controls'][key]['Refine'] = refine
4321
4322    def set_formula(self,*args):
4323        '''Set the chemical formula for the PDF computation.
4324        Use pdf.set_formula(['Si',1],['O',2]) for SiO2.
4325
4326        :param list item1: The element symbol and number of atoms in formula for first element
4327        :param list item2: The element symbol and number of atoms in formula for second element,...
4328
4329        repeat parameters as needed for all elements in the formula.
4330        '''
4331        powderHist = self.proj.histogram(self.data['PDF Controls']['Sample']['Name'])
4332        inst = powderHist.data['Instrument Parameters'][0]
4333        ElList = self.data['PDF Controls']['ElList']
4334        ElList.clear()
4335        sumVol = 0.
4336        for elem,mult in args:
4337            ElList[elem] = G2elem.GetElInfo(elem,inst)
4338            ElList[elem]['FormulaNo'] = mult
4339            Avol = (4.*np.pi/3.)*ElList[elem]['Drad']**3
4340            sumVol += Avol*ElList[elem]['FormulaNo']
4341        self.data['PDF Controls']['Form Vol'] = max(10.0,sumVol)
4342       
4343    def calculate(self,xydata=None,limits=None,inst=None):
4344        '''Compute the PDF using the current parameters. Results are set
4345        in the PDF object arrays (self.data['PDF Controls']['G(R)'] etc.).
4346        Note that if ``xydata``, is specified, the background histograms(s)
4347        will not be accessed from the project file associated with the current
4348        PDF entry. If ``limits`` and ``inst`` are both specified, no histograms
4349        need be in the current project. However, the self.data['PDF Controls']
4350        sections ('Sample', 'Sample Bkg.','Container Bkg.') must be 
4351        non-blank for the corresponding items to be used from``xydata``.
4352
4353        :param dict xydata: an array containing the Sample's I vs Q, and
4354          any or none of the Sample Background, the Container scattering and
4355          the Container Background. If xydata is None (default), the values are
4356          taken from histograms, as named in the PDF's self.data['PDF Controls']
4357          entries with keys 'Sample', 'Sample Bkg.','Container Bkg.' &
4358          'Container'.
4359        :param list limits: upper and lower Q values to be used for PDF
4360          computation. If None (default), the values are
4361          taken from the Sample histogram's .data['Limits'][1] values.
4362        :param dict inst: The Sample histogram's instrument parameters
4363          to be used for PDF computation. If None (default), the values are
4364          taken from the Sample histogram's .data['Instrument Parameters'][0]
4365          values.
4366        '''
4367        data = self.data['PDF Controls']
4368        if xydata is None:
4369            xydata = {}
4370            for key in 'Sample Bkg.','Container Bkg.','Container','Sample':
4371                name = data[key]['Name'].strip()
4372                if name:
4373                    xydata[key] = self.proj.histogram(name).data['data']
4374        if limits is None:
4375            name = data['Sample']['Name'].strip()
4376            limits = self.proj.histogram(name).data['Limits'][1]
4377        if inst is None:
4378            name = data['Sample']['Name'].strip()
4379            inst = self.proj.histogram(name).data['Instrument Parameters'][0]
4380        G2pwd.CalcPDF(data,inst,limits,xydata)
4381        data['I(Q)'] = xydata['IofQ']
4382        data['S(Q)'] = xydata['SofQ']
4383        data['F(Q)'] = xydata['FofQ']
4384        data['G(R)'] = xydata['GofR']
4385
4386    def optimize(self,showFit=True,maxCycles=5,
4387                     xydata=None,limits=None,inst=None):
4388        '''Optimize the low R portion of G(R) to minimize selected
4389        parameters. Note that this updates the parameters in the settings
4390        (self.data['PDF Controls']) but does not update the PDF object
4391        arrays (self.data['PDF Controls']['G(R)'] etc.) with the computed
4392        values, use :meth:`calculate` after a fit to do that.
4393
4394        :param bool showFit: if True (default) the optimized parameters
4395          are shown before and after the fit, as well as the RMS value
4396          in the minimized region.
4397        :param int maxCycles: the maximum number of least-squares cycles;
4398          defaults to 5.
4399        :returns: the result from the optimizer as True or False, depending
4400          on if the refinement converged.
4401        :param dict xydata: an array containing the Sample's I vs Q, and
4402          any or none of the Sample Background, the Container scattering and
4403          the Container Background. If xydata is None (default), the values are
4404          taken from histograms, as named in the PDF's self.data['PDF Controls']
4405          entries with keys 'Sample', 'Sample Bkg.','Container Bkg.' &
4406          'Container'.
4407        :param list limits: upper and lower Q values to be used for PDF
4408          computation. If None (default), the values are
4409          taken from the Sample histogram's .data['Limits'][1] values.
4410        :param dict inst: The Sample histogram's instrument parameters
4411          to be used for PDF computation. If None (default), the values are
4412          taken from the Sample histogram's .data['Instrument Parameters'][0]
4413          values.
4414        '''
4415        data = self.data['PDF Controls']
4416        if xydata is None:
4417            xydata = {}
4418            for key in 'Sample Bkg.','Container Bkg.','Container','Sample':
4419                name = data[key]['Name'].strip()
4420                if name:
4421                    xydata[key] = self.proj.histogram(name).data['data']
4422        if limits is None:
4423            name = data['Sample']['Name'].strip()
4424            limits = self.proj.histogram(name).data['Limits'][1]
4425        if inst is None:
4426            name = data['Sample']['Name'].strip()
4427            inst = self.proj.histogram(name).data['Instrument Parameters'][0]
4428
4429        res = G2pwd.OptimizePDF(data,xydata,limits,inst,showFit,maxCycles)
4430        return res['success']
4431       
4432    def export(self,fileroot,formats):
4433        '''Write out the PDF-related data (G(r), S(Q),...) into files
4434
4435        :param str fileroot: name of file(s) to be written. The extension
4436          will be ignored and set to .iq, .sq, .fq or .gr depending
4437          on the formats selected.
4438        :param str formats: string specifying the file format(s) to be written,
4439          should contain at least one of the following keywords:
4440          I(Q), S(Q), F(Q), G(r) and/or PDFgui (capitalization and
4441          punctuation is ignored). Note that G(r) and PDFgui should not
4442          be specifed together.
4443        '''
4444        PDFsaves = 5*[False]
4445        PDFentry = self.name
4446        name = self.data['PDF Controls']['Sample']['Name'].strip()
4447        limits = self.proj.histogram(name).data['Limits']
4448        inst = self.proj.histogram(name).data['Instrument Parameters'][0]
4449        for i,lbl in enumerate(['I(Q)', 'S(Q)', 'F(Q)', 'G(r)', 'PDFgui']):
4450            PDFsaves[i] = lbl.lower() in formats.lower()
4451        G2fil.PDFWrite(PDFentry,fileroot,PDFsaves,self.data['PDF Controls'],inst,limits)
4452
4453class G2Image(G2ObjectWrapper):
4454    """Wrapper for an IMG tree entry, containing an image and various metadata.
4455    Note that in a GSASIIscriptable script, instances of G2Image will be created by
4456    calls to :meth:`G2Project.add_image` or :meth:`G2Project.images`, not via calls
4457    to :meth:`G2Image.__init__`.
4458
4459    Example use of G2Image:
4460
4461    >>> gpx = G2sc.G2Project(filename='itest.gpx')
4462    >>> imlst = gpx.add_image(idata,fmthint="TIF")
4463    >>> imlst[0].loadControls('stdSettings.imctrl')
4464    >>> imlst[0].setCalibrant('Si    SRM640c')
4465    >>> imlst[0].loadMasks('stdMasks.immask')
4466    >>> imlst[0].Recalibrate()
4467    >>> imlst[0].setControl('outAzimuths',3)
4468    >>> pwdrList = imlst[0].Integrate()
4469
4470
4471    Sample script using an image:   
4472
4473.. code-block::  python
4474
4475    import os,sys
4476    sys.path.insert(0,'/Users/toby/software/G2/GSASII')
4477    import GSASIIscriptable as G2sc
4478    datadir = "/tmp/V4343_IntegrationError/"
4479    PathWrap = lambda fil: os.path.join(datadir,fil)
4480
4481    gpx = G2sc.G2Project(filename=PathWrap('inttest.gpx'))
4482    imlst = gpx.add_image(PathWrap('Si_free_dc800_1-00000.tif'),fmthint="TIF")
4483    imlst[0].loadControls(PathWrap('Si_free_dc800_1-00000.imctrl'))
4484    pwdrList = imlst[0].Integrate()
4485    gpx.save()
4486
4487    """
4488    # parameters in that can be accessed via setControl. This may need future attention
4489    ControlList = {
4490        'int': ['calibskip', 'pixLimit', 'edgemin', 'outChannels',
4491                    'outAzimuths'],
4492        'float': ['cutoff', 'setdist', 'wavelength', 'Flat Bkg',
4493                      'azmthOff', 'tilt', 'calibdmin', 'rotation',
4494                      'distance', 'DetDepth'],
4495        'bool': ['setRings', 'setDefault', 'centerAzm', 'fullIntegrate',
4496                     'DetDepthRef', 'showLines'],
4497        'str': ['SampleShape', 'binType', 'formatName', 'color',
4498                    'type', ],
4499        'list': ['GonioAngles', 'IOtth', 'LRazimuth', 'Oblique', 'PolaVal',
4500                   'SampleAbs', 'center', 'ellipses', 'linescan',
4501                    'pixelSize', 'range', 'ring', 'rings', 'size', ],
4502        'dict': ['varyList'],
4503        }
4504    '''Defines the items known to exist in the Image Controls tree section
4505    and the item's data types. A few are not included here
4506    ('background image', 'dark image', 'Gain map', and 'calibrant') because
4507    these items have special set routines,
4508    where references to entries are checked to make sure their values are
4509    correct.
4510    ''' 
4511       
4512    def __init__(self, data, name, proj):
4513        self.data = data
4514        self.name = name
4515        self.proj = proj
4516
4517    def setControl(self,arg,value):
4518        '''Set an Image Controls parameter in the current image.
4519        If the parameter is not found an exception is raised.
4520
4521        :param str arg: the name of a parameter (dict entry) in the
4522          image. The parameter must be found in :data:`ControlList`
4523          or an exception is raised.
4524        :param value: the value to set the parameter. The value is
4525          cast as the appropriate type from :data:`ControlList`.
4526        '''
4527        for typ in self.ControlList:
4528            if arg in self.ControlList[typ]: break
4529        else:
4530            G2fil.G2Print('Allowed args:\n',[nam for nam,typ in self.findControl('')])
4531            raise Exception('arg {} not defined in G2Image.setControl'
4532                                .format(arg))
4533        try:
4534            if typ == 'int':
4535                self.data['Image Controls'][arg] = int(value)
4536            elif typ == 'float':
4537                self.data['Image Controls'][arg] = float(value)
4538            elif typ == 'bool':
4539                self.data['Image Controls'][arg] = bool(value)
4540            elif typ == 'str':
4541                self.data['Image Controls'][arg] = str(value)
4542            elif typ == 'list':
4543                self.data['Image Controls'][arg] = list(value)
4544            elif typ == 'dict':
4545                self.data['Image Controls'][arg] = dict(value)
4546            else:
4547                raise Exception('Unknown type {} for arg {} in  G2Image.setControl'
4548                                    .format(typ,arg))
4549        except:
4550            raise Exception('Error formatting value {} as type {} for arg {} in  G2Image.setControl'
4551                                    .format(value,typ,arg))
4552
4553    def getControl(self,arg):
4554        '''Return an Image Controls parameter in the current image.
4555        If the parameter is not found an exception is raised.
4556
4557        :param str arg: the name of a parameter (dict entry) in the
4558          image.
4559        :returns: the value as a int, float, list,...
4560        '''
4561        if arg in self.data['Image Controls']:
4562            return self.data['Image Controls'][arg]
4563        G2fil.G2Print(self.findControl(''))
4564        raise Exception('arg {} not defined in G2Image.getControl'.format(arg))
4565
4566    def findControl(self,arg=''):
4567        '''Finds the Image Controls parameter(s) in the current image
4568        that match the string in arg. Default is '' which returns all
4569        parameters.
4570
4571            Example:
4572
4573            >>> findControl('calib')
4574            [['calibskip', 'int'], ['calibdmin', 'float'], ['calibrant', 'str']]
4575
4576        :param str arg: a string containing part of the name of a
4577          parameter (dict entry) in the image's Image Controls.
4578        :returns: a list of matching entries in form
4579          [['item','type'], ['item','type'],...] where each 'item' string
4580          contains the sting in arg.
4581        '''
4582        matchList = []
4583        for typ in self.ControlList:
4584            for item in self.ControlList[typ]:
4585                if arg in item:
4586                    matchList.append([item,typ])
4587        return matchList
4588
4589    def setCalibrant(self,calib):
4590        '''Set a calibrant for the current image
4591
4592        :param str calib: specifies a calibrant name which must be one of
4593          the entries in file ImageCalibrants.py. This is validated and
4594          an error provides a list of valid choices.
4595        '''
4596        import ImageCalibrants as calFile
4597        if calib in calFile.Calibrants.keys():
4598            self.data['Image Controls']['calibrant'] = calib
4599            return
4600        G2fil.G2Print('Calibrant {} is not valid. Valid calibrants'.format(calib))
4601        for i in calFile.Calibrants.keys():
4602            if i: G2fil.G2Print('\t"{}"'.format(i))
4603       
4604    def setControlFile(self,typ,imageRef,mult=None):
4605        '''Set a image to be used as a background/dark/gain map image
4606
4607        :param str typ: specifies image type, which must be one of:
4608           'background image', 'dark image', 'gain map'; N.B. only the first
4609           four characters must be specified and case is ignored.
4610        :param imageRef: A reference to the desired image. Either the Image
4611          tree name (str), the image's index (int) or
4612          a image object (:class:`G2Image`)
4613        :param float mult: a multiplier to be applied to the image (not used
4614          for 'Gain map'; required for 'background image', 'dark image'
4615        '''
4616        if 'back' in typ.lower():
4617            key = 'background image'
4618            mult = float(mult)
4619        elif 'dark' in typ.lower():
4620            key = 'dark image'
4621            mult = float(mult)
4622        elif 'gain' in typ.lower():
4623            #key = 'Gain map'
4624            if mult is not None:
4625                G2fil.G2Print('Warning: Ignoring multiplier for Gain map')
4626            mult = None
4627        else:
4628            raise Exception("Invalid typ {} for setControlFile".format(typ))
4629        imgNam = self.proj.image(imageRef).name
4630        if mult is None:
4631            self.data['Image Controls']['Gain map'] = imgNam
4632        else:
4633            self.data['Image Controls'][key] = [imgNam,mult]
4634
4635    def loadControls(self,filename=None,imgDict=None):
4636        '''load controls from a .imctrl file
4637
4638        :param str filename: specifies a file to be read, which should end
4639          with .imctrl (defaults to None, meaning parameters are input
4640          with imgDict.)
4641        :param dict imgDict: contains a set of image parameters (defaults to
4642          None, meaning parameters are input with filename.)
4643        '''
4644        if filename:
4645            File = open(filename,'r')
4646            Slines = File.readlines()
4647            File.close()
4648            G2fil.LoadControls(Slines,self.data['Image Controls'])
4649            G2fil.G2Print('file {} read into {}'.format(filename,self.name))
4650        elif imgDict:
4651            self.data['Image Controls'].update(imgDict)
4652            G2fil.G2Print('Image controls set')
4653        else:
4654            raise Exception("loadControls called without imgDict or filename specified")
4655
4656    def saveControls(self,filename):
4657        '''write current controls values to a .imctrl file
4658
4659        :param str filename: specifies a file to write, which should end
4660          with .imctrl
4661        '''
4662        G2fil.WriteControls(filename,self.data['Image Controls'])
4663        G2fil.G2Print('file {} written from {}'.format(filename,self.name))
4664
4665    def getControls(self,clean=False):
4666        '''returns current Image Controls as a dict
4667
4668        :param bool clean: causes the calbration information to be deleted
4669        '''
4670        ImageControls = copy.deepcopy(self.data['Image Controls'])
4671        if clean:
4672            ImageControls['showLines'] = True
4673            ImageControls['ring'] = []
4674            ImageControls['rings'] = []
4675            ImageControls['ellipses'] = []
4676            ImageControls['setDefault'] = False
4677            for i in 'range','size','GonioAngles':
4678                if i in ImageControls: del ImageControls[i]
4679        return ImageControls
4680   
4681    def setControls(self,controlsDict):
4682        '''uses dict from :meth:`getControls` to set Image Controls for current image
4683        '''
4684        self.data['Image Controls'].update(copy.deepcopy(controlsDict))
4685   
4686   
4687    def loadMasks(self,filename,ignoreThreshold=False):
4688        '''load masks from a .immask file
4689
4690        :param str filename: specifies a file to be read, which should end
4691          with .immask
4692        :param bool ignoreThreshold: If True, masks are loaded with
4693          threshold masks. Default is False which means any Thresholds
4694          in the file are ignored.
4695        '''
4696        G2fil.readMasks(filename,self.data['Masks'],ignoreThreshold)
4697        G2fil.G2Print('file {} read into {}'.format(filename,self.name))
4698
4699    def initMasks(self):
4700        '''Initialize Masks, including resetting the Thresholds values
4701        '''
4702        self.data['Masks'] = {'Points':[],'Rings':[],'Arcs':[],'Polygons':[],'Frames':[]}
4703        ImageZ = _getCorrImage(Readers['Image'],self.proj,self)
4704        Imin = max(0.,np.min(ImageZ))
4705        Imax = np.max(ImageZ)
4706        self.data['Masks']['Thresholds'] = [(0,Imax),[Imin,Imax]]
4707       
4708    def getMasks(self):
4709        '''load masks from an IMG tree entry
4710        '''
4711        return self.data['Masks']
4712       
4713    def setMasks(self,maskDict,resetThresholds=False):
4714        '''load masks dict (from :meth:`getMasks`) into current IMG record
4715
4716        :param dict maskDict: specifies a dict with image parameters
4717        :param bool resetThresholds: If True, Threshold Masks in the
4718          dict are ignored. The default is False which means Threshold
4719          Masks are retained.
4720        '''
4721        self.data['Masks'] = copy.deepcopy(maskDict)
4722        if resetThresholds:
4723            ImageZ = _getCorrImage(Readers['Image'],self.proj,self)
4724            Imin = max(0.,np.min(ImageZ))
4725            Imax = np.max(ImageZ)
4726            self.data['Masks']['Thresholds'] [(0,Imax),[Imin,Imax]]
4727       
4728    def getVary(self,*args):
4729        '''Return the refinement flag(s) for Image Controls parameter(s)
4730        in the current image.
4731        If the parameter is not found, an exception is raised.
4732
4733        :param str arg: the name of a refinement parameter in the
4734          varyList for the image. The name should be one of
4735          'dep', 'det-X', 'det-Y', 'dist', 'phi', 'tilt', or 'wave'
4736        :param str arg1: the name of a parameter (dict entry) as before,
4737          optional
4738
4739
4740        :returns: a list of bool value(s)
4741        '''
4742        res = []
4743        for arg in args:
4744            if arg in self.data['Image Controls']['varyList']:
4745                res.append(self.data['Image Controls']['varyList'][arg])
4746            else:
4747                raise Exception('arg {} not defined in G2Image.getVary'.format(arg))
4748        return res
4749   
4750    def setVary(self,arg,value):
4751        '''Set a refinement flag for Image Controls parameter in the
4752        current image.
4753        If the parameter is not found an exception is raised.
4754
4755        :param str arg: the name of a refinement parameter in the
4756          varyList for the image. The name should be one of
4757          'dep', 'det-X', 'det-Y', 'dist', 'phi', 'tilt', or 'wave'
4758        :param str arg: the name of a parameter (dict entry) in the
4759          image. The parameter must be found in :data:`ControlList`
4760          or an exception is raised.
4761        :param value: the value to set the parameter. The value is
4762          cast as the appropriate type from :data:`ControlList`.
4763        '''
4764        if arg in self.data['Image Controls']['varyList']:
4765            self.data['Image Controls']['varyList'][arg] = bool(value)
4766        else:
4767            raise Exception('arg {} not defined in G2Image.setVary'.format(arg))
4768
4769    def Recalibrate(self):
4770        '''Invokes a recalibration fit (same as Image Controls/Calibration/Recalibrate
4771        menu command). Note that for this to work properly, the calibration
4772        coefficients (center, wavelength, ditance & tilts) must be fairly close.
4773        This may produce a better result if run more than once.
4774        '''
4775        LoadG2fil()
4776        ImageZ = _getCorrImage(Readers['Image'],self.proj,self)
4777        G2img.ImageRecalibrate(None,ImageZ,self.data['Image Controls'],self.data['Masks'])
4778
4779    def Integrate(self,name=None):
4780        '''Invokes an image integration (same as Image Controls/Integration/Integrate
4781        menu command). All parameters will have previously been set with Image Controls
4782        so no input is needed here. Note that if integration is performed on an
4783        image more than once, histogram entries may be overwritten. Use the name
4784        parameter to prevent this if desired.
4785
4786        :param str name: base name for created histogram(s). If None (default),
4787          the histogram name is taken from the image name.
4788        :returns: a list of created histogram (:class:`G2PwdrData`) objects.
4789        '''
4790        blkSize = 256   #256 seems to be optimal; will break in polymask if >1024
4791        ImageZ = _getCorrImage(Readers['Image'],self.proj,self)
4792        # do integration
4793        ints,azms,Xvals,cancel = G2img.ImageIntegrate(ImageZ,self.data['Image Controls'],self.data['Masks'],blkSize=blkSize)
4794        # code from here on based on G2IO.SaveIntegration, but places results in the current
4795        # project rather than tree
4796        X = Xvals[:-1]
4797        N = len(X)
4798
4799        data = self.data['Image Controls']
4800        Comments = self.data['Comments']
4801        # make name from image, unless overridden
4802        if name:
4803            if not name.startswith(data['type']+' '):
4804                name = data['type']+' '+name
4805        else:
4806            name = self.name.replace('IMG ',data['type']+' ')
4807        if 'PWDR' in name:
4808            if 'target' in data:
4809                names = ['Type','Lam1','Lam2','I(L2)/I(L1)','Zero','Polariz.','U','V','W','X','Y','Z','SH/L','Azimuth'] 
4810                codes = [0 for i in range(14)]
4811            else:
4812                names = ['Type','Lam','Zero','Polariz.','U','V','W','X','Y','Z','SH/L','Azimuth'] 
4813                codes = [0 for i in range(12)]
4814        elif 'SASD' in name:
4815            names = ['Type','Lam','Zero','Azimuth'] 
4816            codes = [0 for i in range(4)]
4817            X = 4.*np.pi*npsind(X/2.)/data['wavelength']    #convert to q
4818        Xminmax = [X[0],X[-1]]
4819        Azms = []
4820        dazm = 0.
4821        if data['fullIntegrate'] and data['outAzimuths'] == 1:
4822            Azms = [45.0,]                              #a poor man's average?
4823        else:
4824            for i,azm in enumerate(azms[:-1]):
4825                if azm > 360. and azms[i+1] > 360.:
4826                    Azms.append(G2img.meanAzm(azm%360.,azms[i+1]%360.))
4827                else:   
4828                    Azms.append(G2img.meanAzm(azm,azms[i+1]))
4829            dazm = np.min(np.abs(np.diff(azms)))/2.
4830        # pull out integration results and make histograms for each
4831        IntgOutList = []
4832        for i,azm in enumerate(azms[:-1]):
4833            Aname = name+" Azm= %.2f"%((azm+dazm)%360.)
4834            # MT dict to contain histogram
4835            HistDict = {}
4836            histItems = [Aname]
4837            Sample = G2obj.SetDefaultSample()       #set as Debye-Scherrer
4838            Sample['Gonio. radius'] = data['distance']
4839            Sample['Omega'] = data['GonioAngles'][0]
4840            Sample['Chi'] = data['GonioAngles'][1]
4841            Sample['Phi'] = data['GonioAngles'][2]
4842            Sample['Azimuth'] = (azm+dazm)%360.    #put here as bin center
4843            polariz = 0.99    #set default polarization for synchrotron radiation!
4844            for item in Comments:
4845                if 'polariz' in item:
4846                    try:
4847                        polariz = float(item.split('=')[1])
4848                    except:
4849                        polariz = 0.99
4850                for key in ('Temperature','Pressure','Time','FreePrm1','FreePrm2','FreePrm3','Omega',
4851                    'Chi','Phi'):
4852                    if key.lower() in item.lower():
4853                        try:
4854                            Sample[key] = float(item.split('=')[1])
4855                        except:
4856                            pass
4857                # if 'label_prm' in item.lower():
4858                #     for num in ('1','2','3'):
4859                #         if 'label_prm'+num in item.lower():
4860                #             Controls['FreePrm'+num] = item.split('=')[1].strip()
4861            if 'PWDR' in Aname:
4862                if 'target' in data:    #from lab x-ray 2D imaging data
4863                    waves = {'CuKa':[1.54051,1.54433],'TiKa':[2.74841,2.75207],'CrKa':[2.28962,2.29351],
4864                                 'FeKa':[1.93597,1.93991],'CoKa':[1.78892,1.79278],'MoKa':[0.70926,0.713543],
4865                                 'AgKa':[0.559363,0.563775]}
4866                    wave1,wave2 = waves[data['target']]
4867                    parms = ['PXC',wave1,wave2,0.5,0.0,polariz,290.,-40.,30.,6.,-14.,0.0,0.0001,Azms[i]]
4868                else:
4869                    parms = ['PXC',data['wavelength'],0.0,polariz,1.0,-0.10,0.4,0.30,1.0,0.0,0.0001,Azms[i]]
4870            elif 'SASD' in Aname:
4871                Sample['Trans'] = data['SampleAbs'][0]
4872                parms = ['LXC',data['wavelength'],0.0,Azms[i]]
4873            Y = ints[i]
4874            Ymin = np.min(Y)
4875            Ymax = np.max(Y)
4876            W = np.where(Y>0.,1./Y,1.e-6)                    #probably not true
4877            section = 'Comments'
4878            histItems += [section]
4879            HistDict[section] = Comments
4880            section = 'Limits'
4881            histItems += [section]
4882            HistDict[section] = copy.deepcopy([tuple(Xminmax),Xminmax])
4883            if 'PWDR' in Aname:
4884                section = 'Background'
4885                histItems += [section]
4886                HistDict[section] = [['chebyschev',1,3,1.0,0.0,0.0],
4887                    {'nDebye':0,'debyeTerms':[],'nPeaks':0,'peaksList':[]}]
4888            inst = [dict(zip(names,zip(parms,parms,codes))),{}]
4889            for item in inst[0]:
4890                inst[0][item] = list(inst[0][item])
4891            section = 'Instrument Parameters'
4892            histItems += [section]
4893            HistDict[section] = inst
4894            if 'PWDR' in Aname:
4895                section = 'Sample Parameters'
4896                histItems += [section]
4897                HistDict[section] = Sample
4898                section = 'Peak List'
4899                histItems += [section]
4900                HistDict[section] = {'sigDict':{},'peaks':[]}
4901                section = 'Index Peak List'
4902                histItems += [section]
4903                HistDict[section] = [[],[]]
4904                section = 'Unit Cells List'
4905                histItems += [section]
4906                HistDict[section] = []
4907                section = 'Reflection Lists'
4908                histItems += [section]
4909                HistDict[section] = {}
4910            # elif 'SASD' in Aname:             
4911            #     section = 'Substances'
4912            #     histItems += [section]
4913            #     HistDict[section] = G2pdG.SetDefaultSubstances()  # this needs to be moved
4914            #     section = 'Sample Parameters'
4915            #     histItems += [section]
4916            #     HistDict[section] = Sample
4917            #     section = 'Models'
4918            #     histItems += [section]
4919            #     HistDict[section] = G2pdG.SetDefaultSASDModel() # this needs to be moved
4920            valuesdict = {
4921                'wtFactor':1.0,'Dummy':False,'ranId':ran.randint(0,sys.maxsize),'Offset':[0.0,0.0],'delOffset':0.02*Ymax,
4922                'refOffset':-0.1*Ymax,'refDelt':0.1*Ymax,'Yminmax':[Ymin,Ymax]}
4923            # if Aname is already in the project replace it
4924            for j in self.proj.names:
4925                if j[0] == Aname: 
4926                    G2fil.G2Print('Replacing "{}" in project'.format(Aname))
4927                    break
4928            else:
4929                G2fil.G2Print('Adding "{}" to project'.format(Aname))
4930                self.proj.names.append([Aname]+
4931                        [u'Comments',u'Limits',u'Background',u'Instrument Parameters',
4932                         u'Sample Parameters', u'Peak List', u'Index Peak List',
4933                         u'Unit Cells List', u'Reflection Lists'])
4934            HistDict['data'] = [valuesdict,
4935                    [np.array(X),np.array(Y),np.array(W),np.zeros(N),np.zeros(N),np.zeros(N)]]
4936            self.proj.data[Aname] = HistDict
4937            IntgOutList.append(self.proj.histogram(Aname))
4938        return IntgOutList
4939
4940##########################
4941# Command Line Interface #
4942##########################
4943# Each of these takes an argparse.Namespace object as their argument,
4944# representing the parsed command-line arguments for the relevant subcommand.
4945# The argument specification for each is in the subcommands dictionary (see
4946# below)
4947
4948commandhelp={}
4949commandhelp["create"] = "creates a GSAS-II project, optionally adding histograms and/or phases"
4950def create(args):
4951    """Implements the create command-line subcommand. This creates a GSAS-II project, optionally adding histograms and/or phases::
4952
4953  usage: GSASIIscriptable.py create [-h] [-d HISTOGRAMS [HISTOGRAMS ...]]
4954                                  [-i IPARAMS [IPARAMS ...]]
4955                                  [-p PHASES [PHASES ...]]
4956                                  filename
4957                                 
4958positional arguments::
4959
4960  filename              the project file to create. should end in .gpx
4961
4962optional arguments::
4963
4964  -h, --help            show this help message and exit
4965  -d HISTOGRAMS [HISTOGRAMS ...], --histograms HISTOGRAMS [HISTOGRAMS ...]
4966                        list of datafiles to add as histograms
4967  -i IPARAMS [IPARAMS ...], --iparams IPARAMS [IPARAMS ...]
4968                        instrument parameter file, must be one for every
4969                        histogram
4970  -p PHASES [PHASES ...], --phases PHASES [PHASES ...]
4971                        list of phases to add. phases are automatically
4972                        associated with all histograms given.
4973
4974    """
4975    proj = G2Project(gpxfile=args.filename)
4976
4977    hist_objs = []
4978    if args.histograms:
4979        for h,i in zip(args.histograms,args.iparams):
4980            G2fil.G2Print("Adding histogram from",h,"with instparm ",i)
4981            hist_objs.append(proj.add_powder_histogram(h, i))
4982
4983    if args.phases: 
4984        for p in args.phases:
4985            G2fil.G2Print("Adding phase from",p)
4986            proj.add_phase(p, histograms=hist_objs)
4987        G2fil.G2Print('Linking phase(s) to histogram(s):')
4988        for h in hist_objs:
4989            G2fil.G2Print ('   '+h.name)
4990
4991    proj.save()
4992
4993commandhelp["add"] = "adds histograms and/or phases to GSAS-II project"
4994def add(args):
4995    """Implements the add command-line subcommand. This adds histograms and/or phases to GSAS-II project::
4996
4997  usage: GSASIIscriptable.py add [-h] [-d HISTOGRAMS [HISTOGRAMS ...]]
4998                               [-i IPARAMS [IPARAMS ...]]
4999                               [-hf HISTOGRAMFORMAT] [-p PHASES [PHASES ...]]
5000                               [-pf PHASEFORMAT] [-l HISTLIST [HISTLIST ...]]
5001                               filename
5002
5003
5004positional arguments::
5005
5006  filename              the project file to open. Should end in .gpx
5007
5008optional arguments::
5009
5010  -h, --help            show this help message and exit
5011  -d HISTOGRAMS [HISTOGRAMS ...], --histograms HISTOGRAMS [HISTOGRAMS ...]
5012                        list of datafiles to add as histograms
5013  -i IPARAMS [IPARAMS ...], --iparams IPARAMS [IPARAMS ...]
5014                        instrument parameter file, must be one for every
5015                        histogram
5016  -hf HISTOGRAMFORMAT, --histogramformat HISTOGRAMFORMAT
5017                        format hint for histogram import. Applies to all
5018                        histograms
5019  -p PHASES [PHASES ...], --phases PHASES [PHASES ...]
5020                        list of phases to add. phases are automatically
5021                        associated with all histograms given.
5022  -pf PHASEFORMAT, --phaseformat PHASEFORMAT
5023                        format hint for phase import. Applies to all phases.
5024                        Example: -pf CIF
5025  -l HISTLIST [HISTLIST ...], --histlist HISTLIST [HISTLIST ...]
5026                        list of histgram indices to associate with added
5027                        phases. If not specified, phases are associated with
5028                        all previously loaded histograms. Example: -l 2 3 4
5029   
5030    """
5031    proj = G2Project(args.filename)
5032
5033    if args.histograms:
5034        for h,i in zip(args.histograms,args.iparams):
5035            G2fil.G2Print("Adding histogram from",h,"with instparm ",i)
5036            proj.add_powder_histogram(h, i, fmthint=args.histogramformat)
5037
5038    if args.phases: 
5039        if not args.histlist:
5040            histlist = proj.histograms()
5041        else:
5042            histlist = [proj.histogram(i) for i in args.histlist]
5043
5044        for p in args.phases:
5045            G2fil.G2Print("Adding phase from",p)
5046            proj.add_phase(p, histograms=histlist, fmthint=args.phaseformat)
5047           
5048        if not args.histlist:
5049            G2fil.G2Print('Linking phase(s) to all histogram(s)')
5050        else:
5051            G2fil.G2Print('Linking phase(s) to histogram(s):')
5052            for h in histlist:
5053                G2fil.G2Print('   '+h.name)
5054    proj.save()
5055
5056
5057commandhelp["dump"] = "Shows the contents of a GSAS-II project"
5058def dump(args):
5059    """Implements the dump command-line subcommand, which shows the contents of a GSAS-II project::
5060
5061       usage: GSASIIscriptable.py dump [-h] [-d] [-p] [-r] files [files ...]
5062
5063positional arguments::
5064
5065  files
5066
5067optional arguments::
5068
5069  -h, --help        show this help message and exit
5070  -d, --histograms  list histograms in files, overrides --raw
5071  -p, --phases      list phases in files, overrides --raw
5072  -r, --raw         dump raw file contents, default
5073 
5074    """
5075    if not args.histograms and not args.phases:
5076        args.raw = True
5077    if args.raw:
5078        import IPython.lib.pretty as pretty
5079
5080    for fname in args.files:
5081        if args.raw:
5082            proj, nameList = LoadDictFromProjFile(fname)
5083            print("file:", fname)
5084            print("names:", nameList)
5085            for key, val in proj.items():
5086                print(key, ":")
5087                pretty.pprint(val)
5088        else:
5089            proj = G2Project(fname)
5090            if args.histograms:
5091                hists = proj.histograms()
5092                for h in hists:
5093                    print(fname, "hist", h.id, h.name)
5094            if args.phases:
5095                phase_list = proj.phases()
5096                for p in phase_list:
5097                    print(fname, "phase", p.id, p.name)
5098
5099
5100commandhelp["browse"] = "Load a GSAS-II project and then open a IPython shell to browse it"
5101def IPyBrowse(args):
5102    """Load a .gpx file and then open a IPython shell to browse it::
5103
5104  usage: GSASIIscriptable.py browse [-h] files [files ...]
5105
5106positional arguments::
5107
5108  files       list of files to browse
5109
5110optional arguments::
5111
5112  -h, --help  show this help message and exit
5113
5114    """
5115    for fname in args.files:
5116        proj, nameList = LoadDictFromProjFile(fname)
5117        msg = "\nfname {} loaded into proj (dict) with names in nameList".format(fname)
5118        GSASIIpath.IPyBreak_base(msg)
5119        break
5120
5121
5122commandhelp["refine"] = '''
5123Conducts refinements on GSAS-II projects according to a list of refinement
5124steps in a JSON dict
5125'''
5126def refine(args):
5127    """Implements the refine command-line subcommand:
5128    conducts refinements on GSAS-II projects according to a JSON refinement dict::
5129
5130        usage: GSASIIscriptable.py refine [-h] gpxfile [refinements]
5131
5132positional arguments::
5133
5134  gpxfile      the project file to refine
5135  refinements  json file of refinements to apply. if not present refines file
5136               as-is
5137
5138optional arguments::
5139
5140  -h, --help   show this help message and exit
5141 
5142    """
5143    proj = G2Project(args.gpxfile)
5144    if args.refinements is None:
5145        proj.refine()
5146    else:
5147        import json
5148        with open(args.refinements) as refs:
5149            refs = json.load(refs)
5150        if type(refs) is not dict:
5151            raise G2ScriptException("Error: JSON object must be a dict.")
5152        if "code" in refs:
5153            print("executing code:\n|  ",'\n|  '.join(refs['code']))
5154            exec('\n'.join(refs['code']))
5155        proj.do_refinements(refs['refinements'])
5156
5157
5158commandhelp["export"] = "Export phase as CIF"
5159def export(args):
5160    """Implements the export command-line subcommand: Exports phase as CIF::
5161
5162      usage: GSASIIscriptable.py export [-h] gpxfile phase exportfile
5163
5164positional arguments::
5165
5166  gpxfile     the project file from which to export
5167  phase       identifier of phase to export
5168  exportfile  the .cif file to export to
5169
5170optional arguments::
5171
5172  -h, --help  show this help message and exit
5173
5174    """
5175    proj = G2Project(args.gpxfile)
5176    phase = proj.phase(args.phase)
5177    phase.export_CIF(args.exportfile)
5178
5179
5180def _args_kwargs(*args, **kwargs):
5181    return args, kwargs
5182
5183# A dictionary of the name of each subcommand, and a tuple
5184# of its associated function and a list of its arguments
5185# The arguments are passed directly to the add_argument() method
5186# of an argparse.ArgumentParser
5187
5188subcommands = {"create":
5189               (create, [_args_kwargs('filename',
5190                                      help='the project file to create. should end in .gpx'),
5191
5192                         _args_kwargs('-d', '--histograms',
5193                                      nargs='+',
5194                                      help='list of datafiles to add as histograms'),
5195                                     
5196                         _args_kwargs('-i', '--iparams',
5197                                      nargs='+',
5198                                      help='instrument parameter file, must be one'
5199                                           ' for every histogram'
5200                                      ),
5201
5202                         _args_kwargs('-p', '--phases',
5203                                      nargs='+',
5204                                      help='list of phases to add. phases are '
5205                                           'automatically associated with all '
5206                                           'histograms given.')]),
5207               "add": (add, [_args_kwargs('filename',
5208                                      help='the project file to open. Should end in .gpx'),
5209
5210                         _args_kwargs('-d', '--histograms',
5211                                      nargs='+',
5212                                      help='list of datafiles to add as histograms'),
5213                                     
5214                         _args_kwargs('-i', '--iparams',
5215                                      nargs='+',
5216                                      help='instrument parameter file, must be one'
5217                                           ' for every histogram'
5218                                      ),
5219                                     
5220                         _args_kwargs('-hf', '--histogramformat',
5221                                      help='format hint for histogram import. Applies to all'
5222                                           ' histograms'
5223                                      ),
5224
5225                         _args_kwargs('-p', '--phases',
5226                                      nargs='+',
5227                                      help='list of phases to add. phases are '
5228                                           'automatically associated with all '
5229                                           'histograms given.'),
5230
5231                         _args_kwargs('-pf', '--phaseformat',
5232                                      help='format hint for phase import. Applies to all'
5233                                           ' phases. Example: -pf CIF'
5234                                      ),
5235                                     
5236                         _args_kwargs('-l', '--histlist',
5237                                      nargs='+',
5238                                      help='list of histgram indices to associate with added'
5239                                           ' phases. If not specified, phases are'
5240                                           ' associated with all previously loaded'
5241                                           ' histograms. Example: -l 2 3 4')]),
5242                                           
5243               "dump": (dump, [_args_kwargs('-d', '--histograms',
5244                                     action='store_true',
5245                                     help='list histograms in files, overrides --raw'),
5246
5247                               _args_kwargs('-p', '--phases',
5248                                            action='store_true',
5249                                            help='list phases in files, overrides --raw'),
5250
5251                               _args_kwargs('-r', '--raw',
5252                                      action='store_true', help='dump raw file contents, default'),
5253
5254                               _args_kwargs('files', nargs='+')]),
5255
5256               "refine":
5257               (refine, [_args_kwargs('gpxfile', help='the project file to refine'),
5258                         _args_kwargs('refinements',
5259                                      help='JSON file of refinements to apply. if not present'
5260                                           ' refines file as-is',
5261                                      default=None,
5262                                      nargs='?')]),
5263
5264               "export": (export, [_args_kwargs('gpxfile',
5265                                                help='the project file from which to export'),
5266                                   _args_kwargs('phase', help='identifier of phase to export'),
5267                                   _args_kwargs('exportfile', help='the .cif file to export to')]),
5268               "browse": (IPyBrowse, [_args_kwargs('files', nargs='+',
5269                                                   help='list of files to browse')])}
5270
5271
5272def main():
5273    '''The command-line interface for calling GSASIIscriptable as a shell command,
5274    where it is expected to be called as::
5275
5276       python GSASIIscriptable.py <subcommand> <file.gpx> <options>
5277
5278    The following subcommands are defined:
5279
5280        * create, see :func:`create`
5281        * add, see :func:`add`
5282        * dump, see :func:`dump`
5283        * refine, see :func:`refine`
5284        * export, :func:`export`
5285        * browse, see :func:`IPyBrowse`
5286
5287    .. seealso::
5288        :func:`create`
5289        :func:`add`
5290        :func:`dump`
5291        :func:`refine`
5292        :func:`export`
5293        :func:`IPyBrowse`
5294    '''
5295    parser = argparse.ArgumentParser(description=
5296        "Use of "+os.path.split(__file__)[1]+" Allows GSAS-II actions from command line."
5297        )
5298    subs = parser.add_subparsers()
5299
5300    # Create all of the specified subparsers
5301    for name, (func, args) in subcommands.items():
5302        new_parser = subs.add_parser(name,help=commandhelp.get(name),
5303                                     description='Command "'+name+'" '+commandhelp.get(name))
5304        for listargs, kwds in args:
5305            new_parser.add_argument(*listargs, **kwds)
5306        new_parser.set_defaults(func=func)
5307
5308    # Parse and trigger subcommand
5309    result = parser.parse_args()
5310    result.func(result)
5311
5312if __name__ == '__main__':
5313    #fname='/tmp/corundum-template.gpx'
5314    #prj = G2Project(fname)
5315    main()
Note: See TracBrowser for help on using the repository browser.