source: trunk/GSASIIscriptable.py @ 4039

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

turn off refinement for simulation

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