source: branch/2frame/GSASIIobj.py @ 2987

Last change on this file since 2987 was 2987, checked in by odonnell, 4 years ago

added some HAP documentation to GSASIIobj

  • Property svn:eol-style set to native
  • Property svn:keywords set to Date Author Revision URL Id
File size: 127.5 KB
Line 
1# -*- coding: utf-8 -*-
2#GSASIIobj - data objects for GSAS-II
3########### SVN repository information ###################
4# $Date: 2017-08-09 21:39:05 +0000 (Wed, 09 Aug 2017) $
5# $Author: odonnell $
6# $Revision: 2987 $
7# $URL: branch/2frame/GSASIIobj.py $
8# $Id: GSASIIobj.py 2987 2017-08-09 21:39:05Z odonnell $
9########### SVN repository information ###################
10
11'''
12*GSASIIobj: Data objects*
13=========================
14
15This module defines and/or documents the data structures used in GSAS-II, as well
16as provides misc. support routines.
17
18Constraints Tree Item
19----------------------
20
21.. _Constraints_table:
22
23.. index::
24   single: Constraints object description
25   single: Data object descriptions; Constraints
26
27Constraints are stored in a dict, separated into groups.
28Note that parameter are named in the following pattern,
29p:h:<var>:n, where p is the phase number, h is the histogram number
30<var> is a variable name and n is the parameter number.
31If a parameter does not depend on a histogram or phase or is unnumbered, that
32number is omitted.
33Note that the contents of each dict item is a List where each element in the
34list is a :ref:`constraint definition objects <Constraint_definitions_table>`.
35The constraints in this form are converted in
36:func:`GSASIIstrIO.ProcessConstraints` to the form used in :mod:`GSASIImapvars`
37
38The keys in the Constraints dict are:
39
40.. tabularcolumns:: |l|p{4.5in}|
41
42==========  ====================================================
43  key         explanation
44==========  ====================================================
45Hist        This specifies a list of constraints on
46            histogram-related parameters,
47            which will be of form :h:<var>:n.
48HAP         This specifies a list of constraints on parameters
49            that are defined for every histogram in each phase
50            and are of form p:h:<var>:n.           
51Phase       This specifies a list of constraints on phase
52            parameters,
53            which will be of form p::<var>:n.
54Global      This specifies a list of constraints on parameters
55            that are not tied to a histogram or phase and
56            are of form ::<var>:n
57==========  ====================================================
58
59.. _Constraint_definitions_table:
60
61.. index::
62   single: Constraint definition object description
63   single: Data object descriptions; Constraint Definition
64
65Each constraint is defined as an item in a list. Each constraint is of form::
66
67[[<mult1>, <var1>], [<mult2>, <var2>],..., <fixedval>, <varyflag>, <constype>]
68
69Where the variable pair list item containing two values [<mult>, <var>], where:
70
71  * <mult> is a multiplier for the constraint (float)
72  * <var> a :class:`G2VarObj` object (previously a str variable name of form
73      'p:h:name[:at]')
74
75Note that the last three items in the list play a special role:
76
77 * <fixedval> is the fixed value for a `constant equation` (``constype=c``)
78   constraint or is None. For a `New variable` (``constype=f``) constraint,
79   a variable name can be specified as a str (used for externally
80   generated constraints)
81 * <varyflag> is True or False for `New variable` (``constype=f``) constraints
82   or is None. This will be implemented in the future to indicate if these variables
83   should be refined.
84 * <constype> is one of four letters, 'e', 'c', 'h', 'f' that determines the type of constraint:
85
86    * 'e' defines a set of equivalent variables. Only the first variable is refined (if the
87      appropriate refine flag is set) and and all other equivalent variables in the list
88      are generated from that variable, using the appropriate multipliers.
89    * 'c' defines a constraint equation of form,
90      :math:`m_1 \\times var_1 + m_2 \\times var_2 + ... = c`
91    * 'h' defines a variable to hold (not vary). Any variable on this list is not varied,
92      even if its refinement flag is set. Only one [mult,var] pair is allowed in a hold
93      constraint and the mult value is ignored.
94      This is of particular value when needing to hold one or more variables where a
95      single flag controls a set of variables such as, coordinates,
96      the reciprocal metric tensor or anisotropic displacement parameter.
97    * 'f' defines a new variable (function) according to relationship
98      :math:`newvar = m_1 \\times var_1 + m_2 \\times var_2 + ...`
99
100Covariance Tree Item
101--------------------
102
103.. _Covariance_table:
104
105.. index::
106   single: Covariance description
107   single: Data object descriptions; Covariance
108
109The Covariance tree item has results from the last least-squares run. They
110are stored in a dict with these keys:
111
112.. tabularcolumns:: |l|l|p{4in}|
113
114=============  ===============  ====================================================
115  key            sub-key        explanation
116=============  ===============  ====================================================
117newCellDict    \                dict with lattice parameters computed by
118                                :func:`GSASIIstrMath.GetNewCellParms` (dict)
119title          \                Name of gpx file(?) (str)
120variables      \                Values for all N refined variables
121                                (list of float values, length N,
122                                ordered to match varyList)
123sig            \                Uncertainty values for all N refined variables
124                                (list of float values, length N,
125                                ordered to match varyList)
126varyList       \                List of directly refined variables
127                                (list of str values, length N)
128newAtomDict    \                dict with atom position values computed in
129                                :func:`GSASIIstrMath.ApplyXYZshifts` (dict)
130Rvals          \                R-factors, GOF, Marquardt value for last
131                                refinement cycle (dict)
132\              Nobs             Number of observed data points (int)
133\              Rwp              overall weighted profile R-factor (%, float)
134\              chisq            sum[w*(Iobs-Icalc)**2] for all data
135                                note this is not the reduced chi squared (float)
136\              lamMax           Marquardt value applied to Hessian diagonal
137                                (float)
138\              GOF              The goodness-of-fit, aka square root of
139                                the reduced chi squared. (float)
140covMatrix      \                The (NxN) covVariance matrix (np.array)
141=============  ===============  ====================================================
142
143Phase Tree Items
144----------------
145
146.. _Phase_table:
147
148.. index::
149   single: Phase object description
150   single: Data object descriptions; Phase
151
152Phase information is stored in the GSAS-II data tree as children of the
153Phases item in a dict with keys:
154
155.. tabularcolumns:: |l|l|p{4in}|
156
157==========  ===============  ====================================================
158  key         sub-key        explanation
159==========  ===============  ====================================================
160General         \            Overall information for the phase (dict)
161  \         AtomPtrs         list of four locations to use to pull info
162                             from the atom records (list)
163  \         F000X            x-ray F(000) intensity (float)
164  \         F000N            neutron F(000) intensity (float)
165  \         Mydir            directory of current .gpx file (str)
166  \         MCSA controls    Monte Carlo-Simulated Annealing controls (dict)
167  \         Cell             List with 8 items: cell refinement flag (bool)
168                             a, b, c, (Angstrom, float)
169                             alpha, beta & gamma (degrees, float)
170                             volume (A^3, float)
171  \         Type             'nuclear' or 'macromolecular' for now (str)
172  \         Map              dict of map parameters
173  \         SH Texture       dict of spherical harmonic preferred orientation
174                             parameters
175  \         Isotope          dict of isotopes for each atom type
176  \         Isotopes         dict of scattering lengths for each isotope
177                             combination for each element in phase 
178  \         Name             phase name (str)
179  \         SGData           Space group details as a :ref:`space group (SGData) object <SGData_table>`
180                             as defined in :func:`GSASIIspc.SpcGroup`.
181  \         Pawley neg wt    Restraint value for negative Pawley intensities
182                             (float)
183  \         Flip             dict of Charge flip controls
184  \         Data plot type   data plot type ('Mustrain', 'Size' or
185                             'Preferred orientation') for powder data (str)
186  \         Mass             Mass of unit cell contents in g/mol
187  \         POhkl            March-Dollase preferred orientation direction
188  \         Z                dict of atomic numbers for each atom type
189  \         vdWRadii         dict of van der Waals radii for each atom type
190  \         Color            Colors for atoms (list of (r,b,g) triplets)
191  \         AtomTypes        List of atom types
192  \         AtomMass         List of masses for atoms
193  \         doPawley         Flag for Pawley intensity extraction (bool)
194  \         NoAtoms          Number of atoms per unit cell of each type (dict)
195  \         Pawley dmin      maximum Q (as d-space) to use for Pawley
196                             extraction (float)
197  \         BondRadii        Default radius for each atom used to compute
198                             interatomic distances (list of floats)
199  \         AngleRadii       Default radius for each atom used to compute
200                             interatomic angles (list of floats)
201  \         DisAglCtls       Dict with distance/angle search controls,
202                             which has keys 'Name', 'AtomTypes',
203                             'BondRadii', 'AngleRadii' which are as above
204                             except are possibly edited. Also contains
205                             'Factors', which is a 2 element list with
206                             a multiplier for bond and angle search range
207                             [typically (0.85,0.85)].
208ranId           \            unique random number Id for phase (int)
209pId             \            Phase Id number for current project (int).
210Atoms           \            Atoms in phase as a list of lists. The outer list
211                             is for each atom, the inner list contains varying
212                             items depending on the type of phase, see
213                             the :ref:`Atom Records <Atoms_table>` description.
214                             (list of lists)
215Drawing         \            Display parameters (dict)
216\           ballScale        Size of spheres in ball-and-stick display (float)
217\           bondList         dict with bonds
218\           contourLevel     map contour level in e/A^3 (float)
219\           showABC          Flag to show view point triplet (bool). True=show.
220\           viewDir          cartesian viewing direction (np.array with three
221                             elements)
222\           Zclip            clipping distance in A (float)
223\           backColor        background for plot as and R,G,B triplet
224                             (default = [0, 0, 0], black).
225                             (list with three atoms)
226\           selectedAtoms    List of selected atoms (list of int values)
227\           showRigidBodies  Flag to highlight rigid body placement
228\           sizeH            Size ratio for H atoms (float)
229\           bondRadius       Size of binds in A (float)
230\           atomPtrs         positions of x, type, site sym, ADP flag in Draw Atoms (list)
231\           viewPoint        list of lists. First item in list is [x,y,z]
232                             in fractional coordinates for the center of
233                             the plot. Second item list of previous & current
234                             atom number viewed (may be [0,0])
235\           showHydrogen     Flag to control plotting of H atoms.
236\           unitCellBox      Flag to control display of the unit cell.
237\           ellipseProb      Probability limit for display of thermal
238                             ellipsoids in % (float).
239\           vdwScale         Multiplier of van der Waals radius for
240                             display of vdW spheres.
241\           Atoms            A list of lists with an entry for each atom
242                             that is plotted.
243\           Zstep            Step to de/increase Z-clip (float)
244\           Quaternion       Viewing quaternion (4 element np.array)
245\           radiusFactor     Distance ratio for searching for bonds. ? Bonds
246                             are located that are within r(Ra+Rb) and (Ra+Rb)/r
247                             where Ra and Rb are the atomic radii.
248\           oldxy            previous view point (list with two floats)
249\           cameraPos        Viewing position in A for plot (float)
250\           depthFog         True if use depthFog on plot - set currently as False (bool)
251RBModels        \            Rigid body assignments (note Rigid body definitions
252                             are stored in their own main top-level tree entry.)
253Pawley ref      \            Pawley reflections
254Histograms      \            A dict of dicts. The key for the outer dict is
255                             the histograms tied to this phase. The inner
256                             dict contains the combined phase/histogram
257                             parameters for items such as scale factors,
258                             size and strain parameters. The following are the
259                             keys to the inner dict. (dict)
260\           Babinet          Dictionary with two entries, 'BabA', 'BabU'
261\           Extinction       Extinction parameter (list of float, bool)
262\           HStrain          Hydrostatic strain. List of two lists. The first is
263                             a list of the HStrain parameters (1-3 depending on
264                             unit cell), the second is a list of boolean
265                             refinement parameters (same length)
266\           Histogram        The name of the associated histogram (str)
267\           LeBail           Flag for LeBail extraction (bool)
268\           Mustrain         List of microstrain parameters, in order:
269
270                               0. Type, one of u'generalized', u'isotropic',
271                                  u'uniaxial'
272                               1. Isotropic/uniaxial parameters - list of 3 floats
273                               2. Refinement flags - list of 3 bools
274                               3. Microstrain axis - list of 3 ints, [h, k, l]
275                               4. Generalized mustrain parameters - list of 2-6
276                                  floats, depending on space group
277                               5. Generalized refinement flags - list of bools,
278                                  corresponding to the parameters of (4)
279\           Pref.Ori.        Preferred Orientation. List of eight parameters, in order:
280
281                               0. Type, 'MD' for March-Dollase or 'SH' for
282                                  Spherical Harmonics
283                               1. Value, float
284                               2. Refinement flag, bool
285                               3. Preferred direction, list of ints, [h, k, l]
286                               4. int
287                               5. dict
288                               6. list
289                               7. float
290\           Scale            Phase fraction, list of [float, bool]
291\           Show
292\           Use              bool
293\           newLeBail        Whether to perform a new LeBail extraction
294MCSA            \            Monte-Carlo simulated annealing parameters (dict)
295\           
296==========  ===============  ====================================================
297
298Rigid Body Objects
299------------------
300
301.. _RBData_table:
302
303.. index::
304   single: Rigid Body Data description
305   single: Data object descriptions; Rigid Body Data
306   
307Rigid body descriptions are available for two types of rigid bodies: 'Vector'
308and 'Residue'. Vector rigid bodies are developed by a sequence of translations each
309with a refinable magnitude and Residue rigid bodies are described as Cartesian coordinates
310with defined refinable torsion angles.
311
312.. tabularcolumns:: |l|l|p{4in}|
313
314==========  ===============  ====================================================
315  key         sub-key        explanation
316==========  ===============  ====================================================
317Vector      RBId             vector rigid bodies (dict of dict)
318\           AtInfo           Drad, Color: atom drawing radius & color for each atom type (dict)
319\           RBname           Name assigned by user to rigid body (str)
320\           VectMag          vector magnitudes in A (list)
321\           rbXYZ            Cartesian coordinates for Vector rigid body (list of 3 float)
322\           rbRef            3 assigned reference atom nos. in rigid body for origin
323                             definition, use center of atoms flag (list of 3 int & 1 bool)
324\           VectRef          refinement flags for VectMag values (list of bool)
325\           rbTypes          Atom types for each atom in rigid body (list of str)
326\           rbVect           Cartesian vectors for each translation used to build rigid body (list of lists)
327\           useCount         Number of times rigid body is used in any structure (int)
328Residue     RBId             residue rigid bodies (dict of dict)
329\           AtInfo           Drad, Color: atom drawing radius & color for each atom type(dict)
330\           RBname           Name assigned by user to rigid body (str)
331\           rbXYZ            Cartesian coordinates for Residue rigid body (list of 3 float)
332\           rbTypes          Atom types for each atom in rigid body (list of str)
333\           atNames          Names of each atom in rigid body (e.g. C1,N2...) (list of str)
334\           rbRef            3 assigned reference atom nos. in rigid body for origin
335                             definition, use center of atoms flag (list of 3 int & 1 bool)
336\           rbSeq            Orig,Piv,angle,Riding (list): definition of internal rigid body
337                             torsion; origin atom (int), pivot atom (int), torsion angle (float),
338                             riding atoms (list of int)
339\           SelSeq           [int,int] used by SeqSizer to identify objects
340\           useCount         Number of times rigid body is used in any structure (int)
341RBIds           \            unique Ids generated upon creation of each rigid body (dict)
342\           Vector           Ids for each Vector rigid body (list)
343\           Residue          Ids for each Residue rigid body (list)
344==========  ===============  ====================================================
345
346Space Group Objects
347-------------------
348
349.. _SGData_table:
350
351.. index::
352   single: Space Group Data description
353   single: Data object descriptions; Space Group Data
354
355Space groups are interpreted by :func:`GSASIIspc.SpcGroup`
356and the information is placed in a SGdata object
357which is a dict with these keys:
358
359.. tabularcolumns:: |l|p{4.5in}|
360
361==========  ====================================================
362  key         explanation
363==========  ====================================================
364SpGrp       space group symbol (str)
365Laue        one of the following 14 Laue classes:
366            -1, 2/m, mmm, 4/m, 4/mmm, 3R,
367            3mR, 3, 3m1, 31m, 6/m, 6/mmm, m3, m3m (str)
368SGInv       True if centrosymmetric, False if not (bool)
369SGLatt      Lattice centering type. Will be one of
370            P, A, B, C, I, F, R (str)
371SGUniq      unique axis if monoclinic. Will be
372            a, b, or c for monoclinic space groups.
373            Will be blank for non-monoclinic. (str)
374SGCen       Symmetry cell centering vectors. A (n,3) np.array
375            of centers. Will always have at least one row:
376            ``np.array([[0, 0, 0]])``
377SGOps       symmetry operations as a list of form
378            ``[[M1,T1], [M2,T2],...]``
379            where :math:`M_n` is a 3x3 np.array
380            and :math:`T_n` is a length 3 np.array.
381            Atom coordinates are transformed where the
382            Asymmetric unit coordinates [X is (x,y,z)]
383            are transformed using
384            :math:`X^\prime = M_n*X+T_n`
385SGSys       symmetry unit cell: type one of
386            'triclinic', 'monoclinic', 'orthorhombic',
387            'tetragonal', 'rhombohedral', 'trigonal',
388            'hexagonal', 'cubic' (str)
389SGPolax     Axes for space group polarity. Will be one of
390            '', 'x', 'y', 'x y', 'z', 'x z', 'y z',
391            'xyz'. In the case where axes are arbitrary
392            '111' is used (P 1, and ?).
393==========  ====================================================
394
395.. _SSGData_table:
396
397.. index::
398   single: Superspace Group Data description
399   single: Data object descriptions; Superspace Group Data
400
401Superspace groups [3+1] are interpreted by :func:`GSASIIspc.SSpcGroup`
402and the information is placed in a SSGdata object
403which is a dict with these keys:
404
405.. tabularcolumns:: |l|p{4.5in}|
406
407==========  ====================================================
408  key         explanation
409==========  ====================================================
410SSpGrp      superspace group symbol extension to space group
411            symbol, accidental spaces removed
412SSGCen      4D cell centering vectors [0,0,0,0] at least
413SSGOps      4D symmetry operations as [M,T] so that M*x+T = x'
414==========  ====================================================
415
416
417Atom Records
418------------
419
420.. _Atoms_table:
421
422.. index::
423   single: Atoms record description
424   single: Data object descriptions; Atoms record
425
426
427If ``phasedict`` points to the phase information in the data tree, then
428atoms are contained in a list of atom records (list) in
429``phasedict['Atoms']``. Also needed to read atom information
430are four pointers, ``cx,ct,cs,cia = phasedict['General']['atomPtrs']``,
431which define locations in the atom record, as shown below. Items shown are
432always present; additional ones for macromolecular phases are marked 'mm'
433
434.. tabularcolumns:: |l|p{4.5in}|
435
436==============   ====================================================
437location         explanation
438==============   ====================================================
439ct-4              mm - residue number (str)
440ct-3              mm - residue name (e.g. ALA) (str)
441ct-2              mm - chain label (str)
442ct-1              atom label (str)
443ct                atom type (str)
444ct+1              refinement flags; combination of 'F', 'X', 'U' (str)
445cx,cx+1,cx+2      the x,y and z coordinates (3 floats)
446cs                site symmetry (str)
447cs+1              site multiplicity (int)
448cia               ADP flag: Isotropic ('I') or Anisotropic ('A')
449cia+1             Uiso (float)
450cia+2...cia+7     U11, U22, U33, U12, U13, U23 (6 floats)
451atom[cia+8]       unique atom identifier (int)
452
453==============   ====================================================
454
455Drawing Atom Records
456--------------------
457
458.. _Drawing_atoms_table:
459
460.. index::
461   single: Drawing atoms record description
462   single: Data object descriptions; Drawing atoms record
463
464
465If ``phasedict`` points to the phase information in the data tree, then
466drawing atoms are contained in a list of drawing atom records (list) in
467``phasedict['Drawing']['Atoms']``. Also needed to read atom information
468are four pointers, ``cx,ct,cs,ci = phasedict['Drawing']['AtomPtrs']``,
469which define locations in the atom record, as shown below. Items shown are
470always present; additional ones for macromolecular phases are marked 'mm'
471
472.. tabularcolumns:: |l|p{4.5in}|
473
474==============   ====================================================
475location         explanation
476==============   ====================================================
477ct-4              mm - residue number (str)
478ct-3              mm - residue name (e.g. ALA) (str)
479ct-2              mm - chain label (str)
480ct-1              atom label (str)
481ct                atom type (str)
482cx,cx+1,cx+2      the x,y and z coordinates (3 floats)
483cs-1              Sym Op symbol; sym. op number + unit cell id (e.g. '1,0,-1') (str)
484cs                atom drawing style; e.g. 'balls & sticks' (str)
485cs+1              atom label style (e.g. 'name') (str)
486cs+2              atom color (RBG triplet) (int)
487cs+3              ADP flag: Isotropic ('I') or Anisotropic ('A')
488cs+4              Uiso (float)
489cs+5...cs+11      U11, U22, U33, U12, U13, U23 (6 floats)
490ci                unique atom identifier; matches source atom Id in Atom Records (int)
491==============   ====================================================
492
493Powder Diffraction Tree Items
494-----------------------------
495
496.. _Powder_table:
497
498.. index::
499   single: Powder data object description
500   single: Data object descriptions; Powder Data
501
502Every powder diffraction histogram is stored in the GSAS-II data tree
503with a top-level entry named beginning with the string "PWDR ". The
504diffraction data for that information are directly associated with
505that tree item and there are a series of children to that item. The
506routines :func:`GSASIIdataGUI.GSASII.GetUsedHistogramsAndPhasesfromTree`
507and :func:`GSASIIstrIO.GetUsedHistogramsAndPhases` will
508load this information into a dictionary where the child tree name is
509used as a key, and the information in the main entry is assigned
510a key of ``Data``, as outlined below.
511
512.. tabularcolumns:: |p{1in}|p{1in}|p{4in}|
513
514======================  ===============  ====================================================
515  key                      sub-key        explanation
516======================  ===============  ====================================================
517Comments                      \           Text strings extracted from the original powder
518                                          data header. These cannot be changed by the user;
519                                          it may be empty.
520Limits                       \            A list of two two element lists, as [[Ld,Hd],[L,H]]
521                                          where L and Ld are the current and default lowest
522                                          two-theta value to be used and
523                                          where H and Hd are the current and default highest
524                                          two-theta value to be used.
525Reflection Lists              \           A dict with an entry for each phase in the
526                                          histogram. The contents of each dict item
527                                          is a dict containing reflections, as described in
528                                          the :ref:`Powder Reflections <PowderRefl_table>`
529                                          description.
530Instrument Parameters         \           A list containing two dicts where the possible
531                                          keys in each dict are listed below. The value
532                                          for each item is a list containing three values:
533                                          the initial value, the current value and a
534                                          refinement flag which can have a value of
535                                          True, False or 0 where 0 indicates a value that
536                                          cannot be refined. The first and second
537                                          values are floats unless otherwise noted.
538                                          Items in the first dict are noted as [1]
539\                         Lam             Specifies a wavelength in Angstroms [1]
540\                         Lam1            Specifies the primary wavelength in
541                                          Angstrom, when an alpha1, alpha2
542                                          source is used [1]
543\                         Lam2            Specifies the secondary wavelength in
544                                          Angstrom, when an alpha1, alpha2
545                                          source is used [1]
546                          I(L2)/I(L1)     Ratio of Lam2 to Lam1 [1]           
547\                         Type            Histogram type (str) [1]:
548                                           * 'PXC' for constant wavelength x-ray
549                                           * 'PNC' for constant wavelength neutron
550                                           * 'PNT' for time of flight neutron
551\                         Zero            Two-theta zero correction in *degrees* [1]
552\                         Azimuth         Azimuthal setting angle for data recorded
553                                          with differing setting angles [1]
554\                         U, V, W         Cagliotti profile coefficients
555                                          for Gaussian instrumental broadening, where the
556                                          FWHM goes as
557                                          :math:`U \\tan^2\\theta + V \\tan\\theta + W` [1]
558\                         X, Y            Cauchy (Lorentzian) instrumental broadening
559                                          coefficients [1]
560\                         SH/L            Variant of the Finger-Cox-Jephcoat asymmetric
561                                          peak broadening ratio. Note that this is the
562                                          average between S/L and H/L where S is
563                                          sample height, H is the slit height and
564                                          L is the goniometer diameter. [1]
565\                         Polariz.        Polarization coefficient. [1]
566wtFactor                      \           A weighting factor to increase or decrease
567                                          the leverage of data in the histogram (float).
568                                          A value of 1.0 weights the data with their
569                                          standard uncertainties and a larger value
570                                          increases the weighting of the data (equivalent
571                                          to decreasing the uncertainties).
572Sample Parameters             \           Specifies a dict with parameters that describe how
573                                          the data were collected, as listed
574                                          below. Refinable parameters are a list containing
575                                          a float and a bool, where the second value
576                                          specifies if the value is refined, otherwise
577                                          the value is a float unless otherwise noted.
578\                         Scale           The histogram scale factor (refinable)
579\                         Absorption      The sample absorption coefficient as
580                                          :math:`\\mu r` where r is the radius
581                                          (refinable). Only valid for Debye-Scherrer geometry.
582\                         SurfaceRoughA   Surface roughness parameter A as defined by
583                                          Surotti,J. Appl. Cryst, 5,325-331, 1972.(refinable -
584                                          only valid for Bragg-Brentano geometry)                                         
585\                         SurfaceRoughB   Surface roughness parameter B (refinable -
586                                          only valid for Bragg-Brentano geometry)                                         
587\                         DisplaceX,      Sample displacement from goniometer center
588                          DisplaceY       where Y is along the beam direction and
589                                          X is perpendicular. Units are :math:`\\mu m`
590                                          (refinable).
591\                         Phi, Chi,       Goniometer sample setting angles, in degrees.
592                          Omega
593\                         Gonio. radius   Radius of the diffractometer in mm
594\                         InstrName       A name for the instrument, used in preparing
595                                          a CIF (str).
596\                         Force,          Variables that describe how the measurement
597                          Temperature,    was performed. Not used directly in
598                          Humidity,       any computations.
599                          Pressure,
600                          Voltage
601\                         ranId           The random-number Id for the histogram
602                                          (same value as where top-level key is ranId)
603\                         Type            Type of diffraction data, may be 'Debye-Scherrer'
604                                          or 'Bragg-Brentano' (str).
605\                         Diffuse         not in use?
606hId                           \           The number assigned to the histogram when
607                                          the project is loaded or edited (can change)
608ranId                         \           A random number id for the histogram
609                                          that does not change
610Background                    \           The background is stored as a list with where
611                                          the first item in the list is list and the second
612                                          item is a dict. The list contains the background
613                                          function and its coefficients; the dict contains
614                                          Debye diffuse terms and background peaks.
615                                          (TODO: this needs to be expanded.)
616Data                          \           The data consist of a list of 6 np.arrays
617                                          containing in order:
618
619                                           0. the x-postions (two-theta in degrees),
620                                           1. the intensity values (Yobs),
621                                           2. the weights for each Yobs value
622                                           3. the computed intensity values (Ycalc)
623                                           4. the background values
624                                           5. Yobs-Ycalc
625======================  ===============  ====================================================
626
627Powder Reflection Data Structure
628--------------------------------
629
630.. _PowderRefl_table:
631
632.. index::
633   single: Powder reflection object description
634   single: Data object descriptions; Powder Reflections
635   
636For every phase in a histogram, the ``Reflection Lists`` value is a dict
637one element of which is `'RefList'`, which is a np.array containing
638reflections. The columns in that array are documented below.
639
640==========  ====================================================
641  index         explanation
642==========  ====================================================
643 0,1,2       h,k,l (float)
644 3           multiplicity
645 4           d-space, Angstrom
646 5           pos, two-theta
647 6           sig, Gaussian width
648 7           gam, Lorenzian width
649 8           :math:`F_{obs}^2`
650 9           :math:`F_{calc}^2`
651 10          reflection phase, in degrees
652 11          intensity correction for reflection, this times
653             :math:`F_{obs}^2` or :math:`F_{calc}^2` gives Iobs or Icalc
654==========  ====================================================
655
656Single Crystal Tree Items
657-------------------------
658
659.. _Xtal_table:
660
661.. index::
662   single: Single Crystal data object description
663   single: Data object descriptions; Single crystal data
664
665Every single crystal diffraction histogram is stored in the GSAS-II data tree
666with a top-level entry named beginning with the string "HKLF ". The
667diffraction data for that information are directly associated with
668that tree item and there are a series of children to that item. The
669routines :func:`GSASIIdataGUI.GSASII.GetUsedHistogramsAndPhasesfromTree`
670and :func:`GSASIIstrIO.GetUsedHistogramsAndPhases` will
671load this information into a dictionary where the child tree name is
672used as a key, and the information in the main entry is assigned
673a key of ``Data``, as outlined below.
674
675.. tabularcolumns:: |l|l|p{4in}|
676
677======================  ===============  ====================================================
678  key                      sub-key        explanation
679======================  ===============  ====================================================
680Data                          \           A dict that contains the
681                                          reflection table,
682                                          as described in the
683                                          :ref:`Single Crystal Reflections
684                                          <XtalRefl_table>`
685                                          description.
686
687Instrument Parameters         \           A list containing two dicts where the possible
688                                          keys in each dict are listed below. The value
689                                          for most items is a list containing two values:
690                                          the initial value, the current value.
691                                          The first and second
692                                          values are floats unless otherwise noted.
693\                         Lam             Specifies a wavelength in Angstroms (two floats)
694\                         Type            Histogram type (two str values):
695                                           * 'SXC' for constant wavelength x-ray
696                                           * 'SNC' for constant wavelength neutron
697                                           * 'SNT' for time of flight neutron
698\                         InstrName       A name for the instrument, used in preparing
699                                          a CIF (str).
700
701wtFactor                      \           A weighting factor to increase or decrease
702                                          the leverage of data in the histogram (float).
703                                          A value of 1.0 weights the data with their
704                                          standard uncertainties and a larger value
705                                          increases the weighting of the data (equivalent
706                                          to decreasing the uncertainties).
707
708hId                           \           The number assigned to the histogram when
709                                          the project is loaded or edited (can change)
710ranId                         \           A random number id for the histogram
711                                          that does not change
712======================  ===============  ====================================================
713
714Single Crystal Reflection Data Structure
715----------------------------------------
716
717.. _XtalRefl_table:
718
719.. index::
720   single: Single Crystal reflection object description
721   single: Data object descriptions; Single Crystal Reflections
722   
723For every single crystal a histogram, the ``'Data'`` item contains
724the structure factors as an np.array in item `'RefList'`.
725The columns in that array are documented below.
726
727==========  ====================================================
728  index         explanation
729==========  ====================================================
730 0,1,2       h,k,l (float)
731 3           multiplicity
732 4           d-space, Angstrom
733 5           :math:`F_{obs}^2`
734 6           :math:`\sigma(F_{obs}^2)`
735 7           :math:`F_{calc}^2`
736 8           :math:`F_{obs}^2T`
737 9           :math:`F_{calc}^2T`
738 10          reflection phase, in degrees
739 11          intensity correction for reflection, this times
740             :math:`F_{obs}^2` or :math:`F_{calc}^2`
741             gives Iobs or Icalc
742==========  ====================================================
743
744Image Data Structure
745--------------------
746
747.. _Image_table:
748
749.. index::
750   image: Image data object description
751   image: Image object descriptions
752   
753Every 2-dimensional image is stored in the GSAS-II data tree
754with a top-level entry named beginning with the string "IMG ". The
755image data are directly associated with that tree item and there
756are a series of children to that item. The routines :func:`GSASIIdataGUI.GSASII.GetUsedHistogramsAndPhasesfromTree`
757and :func:`GSASIIstrIO.GetUsedHistogramsAndPhases` will
758load this information into a dictionary where the child tree name is
759used as a key, and the information in the main entry is assigned
760a key of ``Data``, as outlined below.
761
762.. tabularcolumns:: |l|l|p{4in}|
763
764======================  ======================  ====================================================
765  key                      sub-key              explanation
766======================  ======================  ====================================================
767Comments                       \                Text strings extracted from the original image data
768                                                header or a metafile. These cannot be changed by 
769                                                the user; it may be empty.                                               
770Image Controls              azmthOff            (float) The offset to be applied to an azimuthal
771                                                value. Accomodates
772                                                detector orientations other than with the detector
773                                                X-axis
774                                                horizontal.
775\                           background image    (list:str,float) The name of a tree item ("IMG ...") that is to be subtracted
776                                                during image integration multiplied by value. It must have the same size/shape as
777                                                the integrated image. NB: value < 0 for subtraction.
778\                           calibrant           (str) The material used for determining the position/orientation
779                                                of the image. The data is obtained from :func:`ImageCalibrants`
780                                                and UserCalibrants.py (supplied by user).
781\                           calibdmin           (float) The minimum d-spacing used during the last calibration run.
782\                           calibskip           (int) The number of expected diffraction lines skipped during the last
783                                                calibration run.
784\                           center              (list:floats) The [X,Y] point in detector coordinates (mm) where the direct beam
785                                                strikes the detector plane as determined by calibration. This point
786                                                does not have to be within the limits of the detector boundaries.
787\                           centerAzm           (bool) If True then the azimuth reported for the integrated slice
788                                                of the image is at the center line otherwise it is at the leading edge.
789\                           color               (str) The name of the colormap used to display the image. Default = 'Paired'.
790\                           cutoff              (float) The minimum value of I/Ib for a point selected in a diffraction ring for
791                                                calibration calculations. See pixLimit for details as how point is found.           
792\                           DetDepth            (float) Coefficient for penetration correction to distance; accounts for diffraction
793                                                ring offset at higher angles. Optionally determined by calibration.
794\                           DetDepthRef         (bool) If True then refine DetDepth during calibration/recalibration calculation.
795\                           distance            (float) The distance (mm) from sample to detector plane.
796\                           ellipses            (list:lists) Each object in ellipses is a list [center,phi,radii,color] where
797                                                center (list) is location (mm) of the ellipse center on the detector plane, phi is the
798                                                rotation of the ellipse minor axis from the x-axis, and radii are the minor & major
799                                                radii of the ellipse. If radii[0] is negative then parameters describe a hyperbola. Color
800                                                is the selected drawing color (one of 'b', 'g' ,'r') for the ellipse/hyperbola.
801\                           edgemin             (float) Not used;  parameter in EdgeFinder code.
802\                           fullIntegrate       (bool) If True then integrate over full 360 deg azimuthal range.
803\                           GonioAngles         (list:floats) The 'Omega','Chi','Phi' goniometer angles used for this image.
804                                                Required for texture calculations.
805\                           invert_x            (bool) If True display the image with the x-axis inverted.
806\                           invert_y            (bool) If True display the image with the y-axis inverted.
807\                           IOtth               (list:floats) The minimum and maximum 2-theta values to be used for integration.
808\                           LRazimuth           (list:floats) The minimum and maximum azimuth values to be used for integration.
809\                           Oblique             (list:float,bool) If True apply a detector absorption correction using the value to the
810                                                intensities obtained during integration.
811\                           outAzimuths         (int) The number of azimuth pie slices.
812\                           outChannels         (int) The number of 2-theta steps.
813\                           pixelSize           (list:ints) The X,Y dimensions (microns) of each pixel.
814\                           pixLimit            (int) A box in the image with 2*pixLimit+1 edges is searched to find the maximum.
815                                                This value (I) along with the minimum (Ib) in the box is reported by :func:`GSASIIimage.ImageLocalMax`
816                                                and subject to cutoff in :func:`GSASIIimage.makeRing`.
817                                                Locations are used to construct rings of points for calibration calcualtions.
818\                           PolaVal             (list:float,bool) If type='SASD' and if True, apply polarization correction to intensities from
819                                                integration using value.
820\                           rings               (list:lists) Each entry is [X,Y,dsp] where X & Y are lists of x,y coordinates around a
821                                                diffraction ring with the same d-spacing (dsp)
822\                           ring                (list) The x,y coordinates of the >5 points on an inner ring
823                                                selected by the user,
824\                           Range               (list) The minimum & maximum values of the image
825\                           rotation            (float) The angle between the x-axis and the vector about which the
826                                                detector is tilted. Constrained to -180 to 180 deg.     
827\                           SampleShape         (str) Currently only 'Cylinder'. Sample shape for Debye-Scherrer experiments; used for absorption
828                                                calculations.
829\                           SampleAbs           (list: float,bool) Value of absorption coefficient for Debye-Scherrer experimnents, flag if True
830                                                to cause correction to be applied.
831\                           setDefault          (bool) If True the use the image controls values for all new images to be read. (might be removed)
832\                           setRings            (bool) If True then display all the selected x,y ring positions (vida supra rings) used in the calibration.           
833\                           showLines           (bool) If True then isplay the integration limits to be used.
834\                           size                (list:int) The number of pixels on the image x & y axes
835\                           type                (str) One of 'PWDR', 'SASD' or 'REFL' for powder, small angle or reflectometry data, respectively.
836\                           tilt                (float) The angle the detector normal makes with the incident beam; range -90 to 90.
837\                           wavelength          (float) Tha radiation wavelength (Angstroms) as entered by the user (or someday obtained from the image header).
838                                               
839Masks                       Arcs                (list: lists) Each entry [2-theta,[azimuth[0],azimuth[1]],thickness] describes an arc mask
840                                                to be excluded from integration
841\                           Frames              (list:lists) Each entry describes the x,y points (3 or more - mm) that describe a frame outside
842                                                of which is excluded from recalibration and integration. Only one frame is allowed.
843\                           Points              (list:lists) Each entry [x,y,radius] (mm) describes an excluded spot on the image to be excluded
844                                                from integration.
845\                           Polygons            (list:lists) Each entry is a list of 3+ [x,y] points (mm) that describe a polygon on the image
846                                                to be excluded from integration.
847\                           Rings               (list: lists) Each entry [2-theta,thickness] describes a ring mask
848                                                to be excluded from integration.
849\                           Thresholds          (list:[tuple,list]) [(Imin,Imax),[Imin,Imax]] This gives lower and upper limits for points on the image to be included
850                                                in integrsation. The tuple is the image intensity limits and the list are those set by the user.   
851                                               
852Stress/Strain               Sample phi          (float) Sample rotation about vertical axis.
853\                           Sample z            (float) Sample translation from the calibration sample position (for Sample phi = 0)
854                                                These will be restricted by space group symmetry; result of strain fit refinement.
855\                           Type                (str) 'True' or 'Conventional': The strain model used for the calculation.
856\                           d-zero              (list:dict) Each item is for a diffraction ring on the image; all items are from the same phase
857                                                and are used to determine the strain tensor.
858                                                The dictionary items are:
859                                                'Dset': (float) True d-spacing for the diffraction ring; entered by the user.
860                                                'Dcalc': (float) Average calculated d-spacing determined from strain coeff.
861                                                'Emat': (list: float) The strain tensor elements e11, e12 & e22 (e21=e12, rest are 0)
862                                                'Esig': (list: float) Esds for Emat from fitting.
863                                                'pixLimit': (int) Search range to find highest point on ring for each data point
864                                                'cutoff': (float) I/Ib cutoff for searching.
865                                                'ImxyObs': (list: lists) [[X],[Y]] observed points to be used for strain calculations.
866                                                'ImtaObs': (list: lists) [[d],[azm]] transformed via detector calibration from ImxyObs.
867                                                'ImtaCalc': (list: lists [[d],[azm]] calculated d-spacing & azimuth from fit.
868                                               
869======================  ======================  ====================================================
870
871Parameter Dictionary
872-------------------------
873
874.. _parmDict_table:
875
876.. index::
877   single: Parameter dictionary
878
879The parameter dictionary contains all of the variable parameters for the refinement.
880The dictionary keys are the name of the parameter (<phase>:<hist>:<name>:<atom>).
881It is prepared in two ways. When loaded from the tree
882(in :meth:`GSASIIdataGUI.GSASII.MakeLSParmDict` and
883:meth:`GSASIIIO.ExportBaseclass.loadParmDict`),
884the values are lists with two elements: ``[value, refine flag]``
885
886When loaded from the GPX file (in
887:func:`GSASIIstrMain.Refine` and :func:`GSASIIstrMain.SeqRefine`), the value in the
888dict is the actual parameter value (usually a float, but sometimes a
889letter or string flag value (such as I or A for iso/anisotropic).
890
891
892*Classes and routines*
893----------------------
894
895'''
896import re
897import imp
898import random as ran
899import sys
900import os.path as ospath
901import cPickle
902import GSASIIpath
903import GSASIImath as G2mth
904import GSASIIspc as G2spc
905import numpy as np
906
907GSASIIpath.SetVersionNumber("$Revision: 2987 $")
908
909DefaultControls = {
910    'deriv type':'analytic Hessian',
911    'min dM/M':0.001,'shift factor':1.,'max cyc':3,'F**2':False,'SVDtol':1.e-6,
912    'UsrReject':{'minF/sig':0,'MinExt':0.01,'MaxDF/F':100.,'MaxD':500.,'MinD':0.05},
913    'Copy2Next':False,'Reverse Seq':False,'HatomFix':False,
914    'Author':'no name',
915    'FreePrm1':'Sample humidity (%)',
916    'FreePrm2':'Sample voltage (V)',
917    'FreePrm3':'Applied load (MN)',
918    'ShowCell':False,
919    }
920'''Values to be used as defaults for the initial contents of the ``Controls``
921data tree item.
922'''
923def StripUnicode(string,subs='.'):
924    '''Strip non-ASCII characters from strings
925   
926    :param str string: string to strip Unicode characters from
927    :param str subs: character(s) to place into string in place of each
928      Unicode character. Defaults to '.'
929
930    :returns: a new string with only ASCII characters
931    '''
932    s = ''
933    for c in string:
934        if ord(c) < 128:
935            s += c
936        else:
937            s += subs
938    return s.encode('ascii','replace')
939
940def MakeUniqueLabel(lbl,labellist):
941    '''Make sure that every a label is unique against a list by adding
942    digits at the end until it is not found in list.
943
944    :param str lbl: the input label
945    :param list labellist: the labels that have already been encountered
946    :returns: lbl if not found in labellist or lbl with ``_1-9`` (or
947      ``_10-99``, etc.) appended at the end
948    '''
949    lbl = StripUnicode(lbl.strip(),'_')
950    if not lbl: # deal with a blank label
951        lbl = '_1'
952    if lbl not in labellist:
953        labellist.append(lbl)
954        return lbl
955    i = 1
956    prefix = lbl
957    if '_' in lbl:
958        prefix = lbl[:lbl.rfind('_')]
959        suffix = lbl[lbl.rfind('_')+1:]
960        try:
961            i = int(suffix)+1
962        except: # suffix could not be parsed
963            i = 1
964            prefix = lbl
965    while prefix+'_'+str(i) in labellist:
966        i += 1
967    else:
968        lbl = prefix+'_'+str(i)
969        labellist.append(lbl)
970    return lbl
971
972PhaseIdLookup = {}
973'''dict listing phase name and random Id keyed by sequential phase index as a str;
974best to access this using :func:`LookupPhaseName`
975'''
976PhaseRanIdLookup = {}
977'''dict listing phase sequential index keyed by phase random Id;
978best to access this using :func:`LookupPhaseId`
979'''
980HistIdLookup = {}
981'''dict listing histogram name and random Id, keyed by sequential histogram index as a str;
982best to access this using :func:`LookupHistName`
983'''
984HistRanIdLookup = {}
985'''dict listing histogram sequential index keyed by histogram random Id;
986best to access this using :func:`LookupHistId`
987'''
988AtomIdLookup = {}
989'''dict listing for each phase index as a str, the atom label and atom random Id,
990keyed by atom sequential index as a str;
991best to access this using :func:`LookupAtomLabel`
992'''
993AtomRanIdLookup = {}
994'''dict listing for each phase the atom sequential index keyed by atom random Id;
995best to access this using :func:`LookupAtomId`
996'''
997ShortPhaseNames = {}
998'''a dict containing a possibly shortened and when non-unique numbered
999version of the phase name. Keyed by the phase sequential index.
1000'''
1001ShortHistNames = {}
1002'''a dict containing a possibly shortened and when non-unique numbered
1003version of the histogram name. Keyed by the histogram sequential index.
1004'''
1005
1006VarDesc = {}
1007''' This dictionary lists descriptions for GSAS-II variables,
1008as set in :func:`CompileVarDesc`. See that function for a description
1009for how keys and values are written.
1010'''
1011
1012reVarDesc = {}
1013''' This dictionary lists descriptions for GSAS-II variables with
1014the same values as :attr:`VarDesc` except that keys have been compiled as
1015regular expressions. Initialized in :func:`CompileVarDesc`.
1016'''
1017# create a default space group object for P1; N.B. fails when building documentation
1018try: 
1019    P1SGData = G2spc.SpcGroup('P 1')[1] # data structure for default space group
1020except TypeError:
1021    pass
1022
1023def GetPhaseNames(fl):
1024    ''' Returns a list of phase names found under 'Phases' in GSASII gpx file
1025    NB: there is another one of these in GSASIIstrIO.py that uses the gpx filename
1026
1027    :param file fl: opened .gpx file
1028    :return: list of phase names
1029    '''
1030    PhaseNames = []
1031    while True:
1032        try:
1033            data = cPickle.load(fl)
1034        except EOFError:
1035            break
1036        datum = data[0]
1037        if 'Phases' == datum[0]:
1038            for datus in data[1:]:
1039                PhaseNames.append(datus[0])
1040    fl.seek(0)          #reposition file
1041    return PhaseNames
1042
1043def SetNewPhase(Name='New Phase',SGData=None,cell=None,Super=None):
1044    '''Create a new phase dict with default values for various parameters
1045
1046    :param str Name: Name for new Phase
1047
1048    :param dict SGData: space group data from :func:`GSASIIspc:SpcGroup`;
1049      defaults to data for P 1
1050
1051    :param list cell: unit cell parameter list; defaults to
1052      [1.0,1.0,1.0,90.,90,90.,1.]
1053
1054    '''
1055    if SGData is None: SGData = P1SGData
1056    if cell is None: cell=[1.0,1.0,1.0,90.,90,90.,1.]
1057    phaseData = {
1058        'ranId':ran.randint(0,sys.maxint),
1059        'General':{
1060            'Name':Name,
1061            'Type':'nuclear',
1062            'Modulated':False,
1063            'AtomPtrs':[3,1,7,9],
1064            'SGData':SGData,
1065            'Cell':[False,]+cell,
1066            'Pawley dmin':1.0,
1067            'Data plot type':'None',
1068            'SH Texture':{
1069                'Order':0,
1070                'Model':'cylindrical',
1071                'Sample omega':[False,0.0],
1072                'Sample chi':[False,0.0],
1073                'Sample phi':[False,0.0],
1074                'SH Coeff':[False,{}],
1075                'SHShow':False,
1076                'PFhkl':[0,0,1],
1077                'PFxyz':[0,0,1],
1078                'PlotType':'Pole figure',
1079                'Penalty':[['',],0.1,False,1.0]}},
1080        'Atoms':[],
1081        'Drawing':{},
1082        'Histograms':{},
1083        'Pawley ref':[],
1084        'RBModels':{},
1085        }
1086    if Super and Super.get('Use',False):
1087        phaseData['General'].update({'Modulated':True,'Super':True,'SuperSg':Super['ssSymb']})
1088        phaseData['General']['SSGData'] = G2spc.SSpcGroup(SGData,Super['ssSymb'])
1089        phaseData['General']['SuperVec'] = [Super['ModVec'],False,Super['maxH']]
1090
1091    return phaseData
1092               
1093def ReadCIF(URLorFile):
1094    '''Open a CIF, which may be specified as a file name or as a URL using PyCifRW
1095    (from James Hester).
1096    The open routine gets confused with DOS names that begin with a letter and colon
1097    "C:\dir\" so this routine will try to open the passed name as a file and if that
1098    fails, try it as a URL
1099
1100    :param str URLorFile: string containing a URL or a file name. Code will try first
1101      to open it as a file and then as a URL.
1102
1103    :returns: a PyCifRW CIF object.
1104    '''
1105    import CifFile as cif # PyCifRW from James Hester
1106
1107    # alternate approach:
1108    #import urllib
1109    #ciffile = 'file:'+urllib.pathname2url(filename)
1110   
1111    try:
1112        fp = open(URLorFile,'r')
1113        cf = cif.ReadCif(fp)
1114        fp.close()
1115        return cf
1116    except IOError:
1117        return cif.ReadCif(URLorFile)
1118
1119def IndexAllIds(Histograms,Phases):
1120    '''Scan through the used phases & histograms and create an index
1121    to the random numbers of phases, histograms and atoms. While doing this,
1122    confirm that assigned random numbers are unique -- just in case lightning
1123    strikes twice in the same place.
1124
1125    Note: this code assumes that the atom random Id (ranId) is the last
1126    element each atom record.
1127
1128    This is called in three places (only): :func:`GSASIIstrIO.GetUsedHistogramsAndPhases`
1129    (which loads the histograms and phases from a GPX file),
1130    :meth:`~GSASIIdataGUI.GSASII.GetUsedHistogramsAndPhasesfromTree`
1131    (which loads the histograms and phases from the data tree.) and
1132    :meth:`GSASIIconstrGUI.UpdateConstraints`
1133    (which displays & edits the constraints in a GUI)
1134
1135    TODO: do we need a lookup for rigid body variables?
1136    '''
1137    # process phases and atoms
1138    PhaseIdLookup.clear()
1139    PhaseRanIdLookup.clear()   
1140    AtomIdLookup.clear()
1141    AtomRanIdLookup.clear()
1142    ShortPhaseNames.clear()
1143    for ph in Phases:
1144        cx,ct,cs,cia = Phases[ph]['General']['AtomPtrs']
1145        ranId = Phases[ph]['ranId'] 
1146        while ranId in PhaseRanIdLookup:
1147            # Found duplicate random Id! note and reassign
1148            print ("\n\n*** Phase "+str(ph)+" has repeated ranId. Fixing.\n")
1149            Phases[ph]['ranId'] = ranId = ran.randint(0,sys.maxint)
1150        pId = str(Phases[ph]['pId'])
1151        PhaseIdLookup[pId] = (ph,ranId)
1152        PhaseRanIdLookup[ranId] = pId
1153        shortname = ph  #[:10]
1154        while shortname in ShortPhaseNames.values():
1155            shortname = ph[:8] + ' ('+ pId + ')'
1156        ShortPhaseNames[pId] = shortname
1157        AtomIdLookup[pId] = {}
1158        AtomRanIdLookup[pId] = {}
1159        for iatm,at in enumerate(Phases[ph]['Atoms']):
1160            ranId = at[cia+8]
1161            while ranId in AtomRanIdLookup[pId]: # check for dups
1162                print ("\n\n*** Phase "+str(ph)+" atom "+str(iatm)+" has repeated ranId. Fixing.\n")
1163                at[cia+8] = ranId = ran.randint(0,sys.maxint)
1164            AtomRanIdLookup[pId][ranId] = str(iatm)
1165            if Phases[ph]['General']['Type'] == 'macromolecular':
1166                label = '%s_%s_%s_%s'%(at[ct-1],at[ct-3],at[ct-4],at[ct-2])
1167            else:
1168                label = at[ct-1]
1169            AtomIdLookup[pId][str(iatm)] = (label,ranId)
1170    # process histograms
1171    HistIdLookup.clear()
1172    HistRanIdLookup.clear()
1173    ShortHistNames.clear()
1174    for hist in Histograms:
1175        ranId = Histograms[hist]['ranId']
1176        while ranId in HistRanIdLookup:
1177            # Found duplicate random Id! note and reassign
1178            print ("\n\n*** Histogram "+str(hist)+" has repeated ranId. Fixing.\n")
1179            Histograms[hist]['ranId'] = ranId = ran.randint(0,sys.maxint)
1180        hId = str(Histograms[hist]['hId'])
1181        HistIdLookup[hId] = (hist,ranId)
1182        HistRanIdLookup[ranId] = hId
1183        shortname = hist[:15]
1184        while shortname in ShortHistNames.values():
1185            shortname = hist[:11] + ' ('+ hId + ')'
1186        ShortHistNames[hId] = shortname
1187
1188def LookupAtomId(pId,ranId):
1189    '''Get the atom number from a phase and atom random Id
1190
1191    :param int/str pId: the sequential number of the phase
1192    :param int ranId: the random Id assigned to an atom
1193
1194    :returns: the index number of the atom (str)
1195    '''
1196    if not AtomRanIdLookup:
1197        raise Exception,'Error: LookupAtomId called before IndexAllIds was run'
1198    if pId is None or pId == '':
1199        raise KeyError,'Error: phase is invalid (None or blank)'
1200    pId = str(pId)
1201    if pId not in AtomRanIdLookup:
1202        raise KeyError,'Error: LookupAtomId does not have phase '+pId
1203    if ranId not in AtomRanIdLookup[pId]:
1204        raise KeyError,'Error: LookupAtomId, ranId '+str(ranId)+' not in AtomRanIdLookup['+pId+']'
1205    return AtomRanIdLookup[pId][ranId]
1206
1207def LookupAtomLabel(pId,index):
1208    '''Get the atom label from a phase and atom index number
1209
1210    :param int/str pId: the sequential number of the phase
1211    :param int index: the index of the atom in the list of atoms
1212
1213    :returns: the label for the atom (str) and the random Id of the atom (int)
1214    '''
1215    if not AtomIdLookup:
1216        raise Exception,'Error: LookupAtomLabel called before IndexAllIds was run'
1217    if pId is None or pId == '':
1218        raise KeyError,'Error: phase is invalid (None or blank)'
1219    pId = str(pId)
1220    if pId not in AtomIdLookup:
1221        raise KeyError,'Error: LookupAtomLabel does not have phase '+pId
1222    if index not in AtomIdLookup[pId]:
1223        raise KeyError,'Error: LookupAtomLabel, ranId '+str(index)+' not in AtomRanIdLookup['+pId+']'
1224    return AtomIdLookup[pId][index]
1225
1226def LookupPhaseId(ranId):
1227    '''Get the phase number and name from a phase random Id
1228
1229    :param int ranId: the random Id assigned to a phase
1230    :returns: the sequential Id (pId) number for the phase (str)
1231    '''
1232    if not PhaseRanIdLookup:
1233        raise Exception,'Error: LookupPhaseId called before IndexAllIds was run'
1234    if ranId not in PhaseRanIdLookup:
1235        raise KeyError,'Error: LookupPhaseId does not have ranId '+str(ranId)
1236    return PhaseRanIdLookup[ranId]
1237
1238def LookupPhaseName(pId):
1239    '''Get the phase number and name from a phase Id
1240
1241    :param int/str pId: the sequential assigned to a phase
1242    :returns:  (phase,ranId) where phase is the name of the phase (str)
1243      and ranId is the random # id for the phase (int)
1244    '''
1245    if not PhaseIdLookup:
1246        raise Exception,'Error: LookupPhaseName called before IndexAllIds was run'
1247    if pId is None or pId == '':
1248        raise KeyError,'Error: phase is invalid (None or blank)'
1249    pId = str(pId)
1250    if pId not in PhaseIdLookup:
1251        raise KeyError,'Error: LookupPhaseName does not have index '+pId
1252    return PhaseIdLookup[pId]
1253
1254def LookupHistId(ranId):
1255    '''Get the histogram number and name from a histogram random Id
1256
1257    :param int ranId: the random Id assigned to a histogram
1258    :returns: the sequential Id (hId) number for the histogram (str)
1259    '''
1260    if not HistRanIdLookup:
1261        raise Exception,'Error: LookupHistId called before IndexAllIds was run'
1262    if ranId not in HistRanIdLookup:
1263        raise KeyError,'Error: LookupHistId does not have ranId '+str(ranId)
1264    return HistRanIdLookup[ranId]
1265
1266def LookupHistName(hId):
1267    '''Get the histogram number and name from a histogram Id
1268
1269    :param int/str hId: the sequential assigned to a histogram
1270    :returns:  (hist,ranId) where hist is the name of the histogram (str)
1271      and ranId is the random # id for the histogram (int)
1272    '''
1273    if not HistIdLookup:
1274        raise Exception,'Error: LookupHistName called before IndexAllIds was run'
1275    if hId is None or hId == '':
1276        raise KeyError,'Error: histogram is invalid (None or blank)'
1277    hId = str(hId)
1278    if hId not in HistIdLookup:
1279        raise KeyError,'Error: LookupHistName does not have index '+hId
1280    return HistIdLookup[hId]
1281
1282def fmtVarDescr(varname):
1283    '''Return a string with a more complete description for a GSAS-II variable
1284
1285    :param str varname: A full G2 variable name with 2 or 3 or 4
1286       colons (<p>:<h>:name[:<a>] or <p>::RBname:<r>:<t>])
1287       
1288    :returns: a string with the description
1289    '''
1290    s,l = VarDescr(varname)
1291    return s+": "+l
1292
1293def VarDescr(varname):
1294    '''Return two strings with a more complete description for a GSAS-II variable
1295
1296    :param str name: A full G2 variable name with 2 or 3 or 4
1297       colons (<p>:<h>:name[:<a>] or <p>::RBname:<r>:<t>])
1298       
1299    :returns: (loc,meaning) where loc describes what item the variable is mapped
1300      (phase, histogram, etc.) and meaning describes what the variable does.
1301    '''
1302   
1303    # special handling for parameter names without a colons
1304    # for now, assume self-defining
1305    if varname.find(':') == -1:
1306        return "Global",varname
1307       
1308    l = getVarDescr(varname)
1309    if not l:
1310        return ("invalid variable name ("+str(varname)+")!"),""
1311#        return "invalid variable name!",""
1312
1313    if not l[-1]:
1314        l[-1] = "(variable needs a definition! Set it in CompileVarDesc)"
1315
1316    if len(l) == 3:         #SASD variable name!
1317        s = 'component:'+l[1]
1318        return s,l[-1]
1319    s = ""
1320    if l[0] is not None and l[1] is not None: # HAP: keep short
1321        if l[2] == "Scale": # fix up ambigous name
1322            l[5] = "Phase fraction"
1323        if l[0] == '*':
1324            lbl = 'Seq. ref.'
1325        else:
1326            lbl = ShortPhaseNames.get(l[0],'? #'+str(l[0]))
1327        if l[1] == '*':
1328            hlbl = 'Seq. ref.'
1329        else:
1330            hlbl = ShortHistNames.get(l[1],'? #'+str(l[1]))
1331        if hlbl[:4] == 'HKLF':
1332            hlbl = 'Xtl='+hlbl[5:]
1333        elif hlbl[:4] == 'PWDR':
1334            hlbl = 'Pwd='+hlbl[5:]
1335        else:
1336            hlbl = 'Hist='+hlbl
1337        s = "Ph="+str(lbl)+" * "+str(hlbl)
1338    else:
1339        if l[2] == "Scale": # fix up ambigous name: must be scale factor, since not HAP
1340            l[5] = "Scale factor"
1341        if l[2] == 'Back': # background parameters are "special", alas
1342            s = 'Hist='+ShortHistNames.get(l[1],'? #'+str(l[1]))
1343            l[-1] += ' #'+str(l[3])
1344        elif l[4] is not None: # rigid body parameter or modulation parm
1345            lbl = ShortPhaseNames.get(l[0],'phase?')
1346            if 'RB' in l[2]:    #rigid body parm
1347                s = "Res #"+str(l[3])+" body #"+str(l[4])+" in "+str(lbl)
1348            else: #modulation parm
1349                s = 'Atom %s wave %s in %s'%(LookupAtomLabel(l[0],l[3])[0],l[4],lbl)
1350        elif l[3] is not None: # atom parameter,
1351            lbl = ShortPhaseNames.get(l[0],'phase?')
1352            try:
1353                albl = LookupAtomLabel(l[0],l[3])[0]
1354            except KeyError:
1355                albl = 'Atom?'
1356            s = "Atom "+str(albl)+" in "+str(lbl)
1357        elif l[0] == '*':
1358            s = "All phases "
1359        elif l[0] is not None:
1360            lbl = ShortPhaseNames.get(l[0],'phase?')
1361            s = "Phase "+str(lbl)
1362        elif l[1] == '*':
1363            s = 'All hists'
1364        elif l[1] is not None:
1365            hlbl = ShortHistNames.get(l[1],'? #'+str(l[1]))
1366            if hlbl[:4] == 'HKLF':
1367                hlbl = 'Xtl='+hlbl[5:]
1368            elif hlbl[:4] == 'PWDR':
1369                hlbl = 'Pwd='+hlbl[5:]
1370            else:
1371                hlbl = 'Hist='+hlbl
1372            s = str(hlbl)
1373    if not s:
1374        s = 'Global'
1375    return s,l[-1]
1376
1377def getVarDescr(varname):
1378    '''Return a short description for a GSAS-II variable
1379
1380    :param str name: A full G2 variable name with 2 or 3 or 4
1381       colons (<p>:<h>:name[:<a1>][:<a2>])
1382     
1383    :returns: a six element list as [`p`,`h`,`name`,`a1`,`a2`,`description`],
1384      where `p`, `h`, `a1`, `a2` are str values or `None`, for the phase number,
1385      the histogram number and the atom number; `name` will always be
1386      a str; and `description` is str or `None`.
1387      If the variable name is incorrectly formed (for example, wrong
1388      number of colons), `None` is returned instead of a list.
1389    '''
1390    l = varname.split(':')
1391    if len(l) == 2:     #SASD parameter name
1392        return varname,l[0],getDescr(l[1])
1393    if len(l) == 3:
1394        l += [None,None]
1395    elif len(l) == 4:
1396        l += [None]
1397    elif len(l) != 5:
1398        return None
1399    for i in (0,1,3,4):
1400        if l[i] == "":
1401            l[i] = None
1402    l += [getDescr(l[2])]
1403    return l
1404   
1405def CompileVarDesc():
1406    '''Set the values in the variable description lookup table (:attr:`VarDesc`)
1407    into :attr:`reVarDesc`. This is called in :func:`getDescr` so the initialization
1408    is always done before use.
1409
1410    Note that keys may contain regular expressions, where '[xyz]'
1411    matches 'x' 'y' or 'z' (equivalently '[x-z]' describes this as range of values).
1412    '.*' matches any string. For example::
1413
1414    'AUiso':'Atomic isotropic displacement parameter',
1415
1416    will match variable ``'p::AUiso:a'``.
1417    If parentheses are used in the key, the contents of those parentheses can be
1418    used in the value, such as::
1419
1420    'AU([123][123])':'Atomic anisotropic displacement parameter U\\1',
1421
1422    will match ``AU11``, ``AU23``,.. and `U11`, `U23` etc will be displayed
1423    in the value when used.
1424   
1425    '''
1426    if reVarDesc: return # already done
1427    for key,value in {
1428        # derived or other sequential vars
1429        '([abc])$' : 'Lattice parameter, \\1, from Ai and Djk', # N.B. '$' prevents match if any characters follow
1430        u'\u03B1' : u'Lattice parameter, \u03B1, from Ai and Djk',
1431        u'\u03B2' : u'Lattice parameter, \u03B2, from Ai and Djk',
1432        u'\u03B3' : u'Lattice parameter, \u03B3, from Ai and Djk',
1433        # ambiguous, alas:
1434        'Scale' : 'Phase or Histogram scale factor',
1435        # Phase vars (p::<var>)
1436        'A([0-5])' : 'Reciprocal metric tensor component \\1',
1437        '[vV]ol' : 'Unit cell volume', # probably an error that both upper and lower case are used
1438        # Atom vars (p::<var>:a)
1439        'dA([xyz])$' : 'change to atomic coordinate, \\1',
1440        'A([xyz])$' : '\\1 fractional atomic coordinate',
1441        'AUiso':'Atomic isotropic displacement parameter',
1442        'AU([123][123])':'Atomic anisotropic displacement parameter U\\1',
1443        'Afrac': 'Atomic site fraction parameter',
1444        'Amul': 'Atomic site multiplicity value',
1445        'AM([xyz])$' : 'Atomic magnetic moment parameter, \\1',
1446        # Hist & Phase (HAP) vars (p:h:<var>)
1447        'Back': 'Background term',
1448        'BkPkint;(.*)':'Background peak #\\1 intensity',
1449        'BkPkpos;(.*)':'Background peak #\\1 position',
1450        'BkPksig;(.*)':'Background peak #\\1 Gaussian width',
1451        'BkPkgam;(.*)':'Background peak #\\1 Cauchy width',
1452        'Bab([AU])': 'Babinet solvent scattering coef. \\1',
1453        'D([123][123])' : 'Anisotropic strain coef. \\1',
1454        'Extinction' : 'Extinction coef.',
1455        'MD' : 'March-Dollase coef.',
1456        'Mustrain;.*' : 'Microstrain coef.',
1457        'Size;.*' : 'Crystallite size value',
1458        'eA$' : 'Cubic mustrain value',
1459        'Ep$' : 'Primary extinction',
1460        'Es$' : 'Secondary type II extinction',
1461        'Eg$' : 'Secondary type I extinction',
1462        'Flack' : 'Flack parameter',
1463        'TwinFr' : 'Twin fraction',
1464        #Histogram vars (:h:<var>)
1465        'Absorption' : 'Absorption coef.',
1466        'Displace([XY])' : 'Debye-Scherrer sample displacement \\1',
1467        'Lam' : 'Wavelength',
1468        'Polariz\.' : 'Polarization correction',
1469        'SH/L' : 'FCJ peak asymmetry correction',
1470        '([UVW])$' : 'Gaussian instrument broadening \\1',
1471        '([XY])$' : 'Cauchy instrument broadening \\1',
1472        'Zero' : 'Debye-Scherrer zero correction',
1473        'nDebye' : 'Debye model background corr. terms',
1474        'nPeaks' : 'Fixed peak background corr. terms',
1475        'RBV.*' : 'Vector rigid body parameter',
1476        'RBR.*' : 'Residue rigid body parameter',
1477        'RBRO([aijk])' : 'Residue rigid body orientation parameter',
1478        'RBRP([xyz])' : 'Residue rigid body position parameter',
1479        'RBRTr;.*' : 'Residue rigid body torsion parameter',
1480        'RBR([TLS])([123AB][123AB])' : 'Residue rigid body group disp. param.',
1481        'constr([0-9]*)' : 'Parameter from constraint',
1482        # supersymmetry parameters  p::<var>:a:o 'Flen','Fcent'?
1483        'mV([0-2])$' : 'Modulation vector component \\1',
1484        'Fsin'  :   'Sin site fraction modulation',
1485        'Fcos'  :   'Cos site fraction modulation',
1486        'Fzero'  :   'Crenel function offset',      #may go away
1487        'Fwid'   :   'Crenel function width',
1488        'Tmin'   :   'ZigZag/Block min location',
1489        'Tmax'   :   'ZigZag/Block max location',
1490        '([XYZ])max': 'ZigZag/Block max value for \\1',
1491        '([XYZ])sin'  : 'Sin position wave for \\1',
1492        '([XYZ])cos'  : 'Cos position wave for \\1',
1493        'U([123][123])sin$' :  'Sin thermal wave for U\\1',
1494        'U([123][123])cos$' :  'Cos thermal wave for U\\1',
1495        'M([XYZ])sin$' :  'Sin mag. moment wave for \\1',
1496        'M([XYZ])cos$' :  'Cos mag. moment wave for \\1',
1497        # PDF peak parms (l:<var>;l = peak no.)
1498        'PDFpos'  : 'PDF peak position',
1499        'PDFmag'  : 'PDF peak magnitude',
1500        'PDFsig'  : 'PDF peak std. dev.',
1501        # SASD vars (l:<var>;l = component)
1502        'Aspect ratio' : 'Particle aspect ratio',
1503        'Length' : 'Cylinder length',
1504        'Diameter' : 'Cylinder/disk diameter',
1505        'Thickness' : 'Disk thickness',
1506        'Shell thickness' : 'Multiplier to get inner(<1) or outer(>1) sphere radius',
1507        'Dist' : 'Interparticle distance',
1508        'VolFr' : 'Dense scatterer volume fraction',
1509        'epis' : 'Sticky sphere epsilon',
1510        'Sticky' : 'Stickyness',
1511        'Depth' : 'Well depth',
1512        'Width' : 'Well width',
1513        'Volume' : 'Particle volume',
1514        'Radius' : 'Sphere/cylinder/disk radius',
1515        'Mean' : 'Particle mean radius',
1516        'StdDev' : 'Standard deviation in Mean',
1517        'G$': 'Guinier prefactor',
1518        'Rg$': 'Guinier radius of gyration',
1519        'B$': 'Porod prefactor',
1520        'P$': 'Porod power',
1521        'Cutoff': 'Porod cutoff',
1522        'PkInt': 'Bragg peak intensity',
1523        'PkPos': 'Bragg peak position',
1524        'PkSig': 'Bragg peak sigma',
1525        'PkGam': 'Bragg peak gamma',
1526        'e([12][12])' : 'strain tensor e\1',   # strain vars e11, e22, e12
1527        'Dcalc': 'Calc. d-spacing',
1528        'Back$': 'background parameter',
1529        'pos$': 'peak position',
1530        'int$': 'peak intensity',
1531        'WgtFrac':'phase weight fraction',
1532        }.items():
1533        VarDesc[key] = value
1534        reVarDesc[re.compile(key)] = value
1535
1536def getDescr(name):
1537    '''Return a short description for a GSAS-II variable
1538
1539    :param str name: The descriptive part of the variable name without colons (:)
1540     
1541    :returns: a short description or None if not found
1542    '''
1543
1544    CompileVarDesc() # compile the regular expressions, if needed
1545    for key in reVarDesc:
1546        m = key.match(name)
1547        if m:
1548            reVarDesc[key]
1549            return m.expand(reVarDesc[key])
1550    return None
1551
1552def GenWildCard(varlist):
1553    '''Generate wildcard versions of G2 variables. These introduce '*'
1554    for a phase, histogram or atom number (but only for one of these
1555    fields) but only when there is more than one matching variable in the
1556    input variable list. So if the input is this::
1557   
1558      varlist = ['0::AUiso:0', '0::AUiso:1', '1::AUiso:0']
1559
1560    then the output will be this::
1561   
1562       wildList = ['*::AUiso:0', '0::AUiso:*']
1563
1564    :param list varlist: an input list of GSAS-II variable names
1565      (such as 0::AUiso:0)
1566
1567    :returns: wildList, the generated list of wild card variable names.
1568    '''
1569    wild = []
1570    for i in (0,1,3):
1571        currentL = varlist[:]
1572        while currentL:
1573            item1 = currentL.pop(0)
1574            i1splt = item1.split(':')
1575            if i >= len(i1splt): continue
1576            if i1splt[i]:
1577                nextL = []
1578                i1splt[i] = '[0-9]+'
1579                rexp = re.compile(':'.join(i1splt))
1580                matchlist = [item1]
1581                for nxtitem in currentL:
1582                    if rexp.match(nxtitem):
1583                        matchlist += [nxtitem]
1584                    else:
1585                        nextL.append(nxtitem)
1586                if len(matchlist) > 1:
1587                    i1splt[i] = '*'
1588                    wild.append(':'.join(i1splt))
1589                currentL = nextL
1590    return wild
1591
1592def LookupWildCard(varname,varlist):
1593    '''returns a list of variable names from list varname
1594    that match wildcard name in varname
1595   
1596    :param str varname: a G2 variable name containing a wildcard
1597      (such as \*::var)
1598    :param list varlist: the list of all variable names used in
1599      the current project
1600    :returns: a list of matching GSAS-II variables (may be empty) 
1601    '''
1602    rexp = re.compile(varname.replace('*','[0-9]+'))
1603    return sorted([var for var in varlist if rexp.match(var)])
1604
1605
1606def _lookup(dic,key):
1607    '''Lookup a key in a dictionary, where None returns an empty string
1608    but an unmatched key returns a question mark. Used in :class:`G2VarObj`
1609    '''
1610    if key is None:
1611        return ""
1612    elif key == "*":
1613        return "*"
1614    else:
1615        return dic.get(key,'?')
1616
1617class G2VarObj(object):
1618    '''Defines a GSAS-II variable either using the phase/atom/histogram
1619    unique Id numbers or using a character string that specifies
1620    variables by phase/atom/histogram number (which can change).
1621    Note that :func:`LoadID` should be used to (re)load the current Ids
1622    before creating or later using the G2VarObj object.
1623
1624    This can store rigid body variables, but does not translate the residue # and
1625    body # to/from random Ids
1626
1627    A :class:`G2VarObj` object can be created with a single parameter:
1628   
1629    :param str/tuple varname: a single value can be used to create a :class:`G2VarObj`
1630      object. If a string, it must be of form "p:h:var" or "p:h:var:a", where
1631
1632     * p is the phase number (which may be left blank or may be '*' to indicate all phases);
1633     * h is the histogram number (which may be left blank or may be '*' to indicate all histograms);
1634     * a is the atom number (which may be left blank in which case the third colon is omitted).
1635       The atom number can be specified as '*' if a phase number is specified (not as '*').
1636       For rigid body variables, specify a will be a string of form "residue:body#"
1637
1638      Alternately a single tuple of form (Phase,Histogram,VarName,AtomID) can be used, where
1639      Phase, Histogram, and AtomID are None or are ranId values (or one can be '*')
1640      and VarName is a string. Note that if Phase is '*' then the AtomID is an atom number.
1641      For a rigid body variables, AtomID is a string of form "residue:body#".
1642
1643    If four positional arguments are supplied, they are:
1644
1645    :param str/int phasenum: The number for the phase (or None or '*')
1646    :param str/int histnum: The number for the histogram (or None or '*')
1647    :param str varname: a single value can be used to create a :class:`G2VarObj`
1648    :param str/int atomnum: The number for the atom (or None or '*')
1649   
1650    '''
1651    IDdict = {}
1652    IDdict['phases'] = {}
1653    IDdict['hists'] = {}
1654    IDdict['atoms'] = {}
1655    def __init__(self,*args):
1656        self.phase = None
1657        self.histogram = None
1658        self.name = ''
1659        self.atom = None
1660        if len(args) == 1 and (type(args[0]) is list or type(args[0]) is tuple) and len(args[0]) == 4:
1661            # single arg with 4 values
1662            self.phase,self.histogram,self.name,self.atom = args[0]
1663        elif len(args) == 1 and ':' in args[0]:
1664            #parse a string
1665            lst = args[0].split(':')
1666            if lst[0] == '*':
1667                self.phase = '*'
1668                if len(lst) > 3:
1669                    self.atom = lst[3]
1670                self.histogram = HistIdLookup.get(lst[1],[None,None])[1]
1671            elif lst[1] == '*':           
1672                self.histogram = '*'
1673                self.phase = PhaseIdLookup.get(lst[0],[None,None])[1]
1674            else:
1675                self.histogram = HistIdLookup.get(lst[1],[None,None])[1]
1676                self.phase = PhaseIdLookup.get(lst[0],[None,None])[1]
1677                if len(lst) == 4:
1678                    if lst[3] == '*':
1679                        self.atom = '*'
1680                    else:
1681                        self.atom = AtomIdLookup[lst[0]].get(lst[3],[None,None])[1]
1682                elif len(lst) == 5:
1683                    self.atom = lst[3]+":"+lst[4]
1684                elif len(lst) == 3:
1685                    pass
1686                else:
1687                    raise Exception,"Too many colons in var name "+str(args[0])
1688            self.name = lst[2]
1689        elif len(args) == 4:
1690            if args[0] == '*':
1691                self.phase = '*'
1692                self.atom = args[3]
1693            else:
1694                self.phase = PhaseIdLookup.get(str(args[0]),[None,None])[1]
1695                if args[3] == '*':
1696                    self.atom = '*'
1697                elif args[0] is not None:
1698                    self.atom = AtomIdLookup[args[0]].get(str(args[3]),[None,None])[1]
1699            if args[1] == '*':
1700                self.histogram = '*'
1701            else:
1702                self.histogram = HistIdLookup.get(str(args[1]),[None,None])[1]
1703            self.name = args[2]
1704        else:
1705            raise Exception,"Incorrectly called GSAS-II parameter name"
1706
1707        #print "DEBUG: created ",self.phase,self.histogram,self.name,self.atom
1708
1709    def __str__(self):
1710        return self.varname()
1711
1712    def varname(self):
1713        '''Formats the GSAS-II variable name as a "traditional" GSAS-II variable
1714        string (p:h:<var>:a) or (p:h:<var>)
1715
1716        :returns: the variable name as a str
1717        '''
1718        a = ""
1719        if self.phase == "*":
1720            ph = "*"
1721            if self.atom:
1722                a = ":" + str(self.atom)
1723        else:
1724            ph = _lookup(PhaseRanIdLookup,self.phase)
1725            if self.atom == '*':
1726                a = ':*'
1727            elif self.atom:
1728                if ":" in str(self.atom):
1729                    a = ":" + str(self.atom)
1730                elif ph in AtomRanIdLookup:
1731                    a = ":" + AtomRanIdLookup[ph].get(self.atom,'?')
1732                else:
1733                    a = ":?"
1734        if self.histogram == "*":
1735            hist = "*"
1736        else:
1737            hist = _lookup(HistRanIdLookup,self.histogram)
1738        s = (ph + ":" + hist + ":" + str(self.name)) + a
1739        return s
1740   
1741    def __repr__(self):
1742        '''Return the detailed contents of the object
1743        '''
1744        s = "<"
1745        if self.phase == '*':
1746            s += "Phases: all; "
1747            if self.atom is not None:
1748                if ":" in str(self.atom):
1749                    s += "Rigid body" + str(self.atom) + "; "
1750                else:
1751                    s += "Atom #" + str(self.atom) + "; "
1752        elif self.phase is not None:
1753            ph =  _lookup(PhaseRanIdLookup,self.phase)
1754            s += "Phase: rId=" + str(self.phase) + " (#"+ ph + "); "
1755            if self.atom == '*':
1756                s += "Atoms: all; "
1757            elif ":" in str(self.atom):
1758                s += "Rigid body" + str(self.atom) + "; "
1759            elif self.atom is not None:
1760                s += "Atom rId=" + str(self.atom)
1761                if ph in AtomRanIdLookup:
1762                    s += " (#" + AtomRanIdLookup[ph].get(self.atom,'?') + "); "
1763                else:
1764                    s += " (#? -- not found!); "
1765        if self.histogram == '*':
1766            s += "Histograms: all; "
1767        elif self.histogram is not None:
1768            hist = _lookup(HistRanIdLookup,self.histogram)
1769            s += "Histogram: rId=" + str(self.histogram) + " (#"+ hist + "); "
1770        s += 'Variable name="' + str(self.name) + '">'
1771        return s+" ("+self.varname()+")"
1772
1773    def __eq__(self, other):
1774        if type(other) is type(self):
1775            return (self.phase == other.phase and
1776                    self.histogram == other.histogram and
1777                    self.name == other.name and
1778                    self.atom == other.atom)
1779        return False
1780
1781    def _show(self):
1782        'For testing, shows the current lookup table'
1783        print 'phases', self.IDdict['phases']
1784        print 'hists', self.IDdict['hists']
1785        print 'atomDict', self.IDdict['atoms']
1786
1787#==========================================================================
1788def SetDefaultSample():
1789    'Fills in default items for the Sample dictionary for Debye-Scherrer & SASD'
1790    return {
1791        'InstrName':'',
1792        'ranId':ran.randint(0,sys.maxint),
1793        'Scale':[1.0,True],'Type':'Debye-Scherrer','Absorption':[0.0,False],
1794        'DisplaceX':[0.0,False],'DisplaceY':[0.0,False],'Diffuse':[],
1795        'Temperature':300.,'Pressure':0.1,'Time':0.0,
1796        'FreePrm1':0.,'FreePrm2':0.,'FreePrm3':0.,
1797        'Gonio. radius':200.0,
1798        'Omega':0.0,'Chi':0.0,'Phi':0.0,'Azimuth':0.0,
1799#SASD items
1800        'Materials':[{'Name':'vacuum','VolFrac':1.0,},{'Name':'vacuum','VolFrac':0.0,}],
1801        'Thick':1.0,'Contrast':[0.0,0.0],       #contrast & anomalous contrast
1802        'Trans':1.0,                            #measured transmission
1803        'SlitLen':0.0,                          #Slit length - in Q(A-1)
1804        }
1805######################################################################
1806class ImportBaseclass(object):
1807    '''Defines a base class for the reading of input files (diffraction
1808    data, coordinates,...). See :ref:`Writing a Import Routine<Import_routines>`
1809    for an explanation on how to use a subclass of this class.
1810    '''
1811    class ImportException(Exception):
1812        '''Defines an Exception that is used when an import routine hits an expected error,
1813        usually in .Reader.
1814
1815        Good practice is that the Reader should define a value in self.errors that
1816        tells the user some information about what is wrong with their file.         
1817        '''
1818        pass
1819   
1820    UseReader = True  # in __init__ set value of self.UseReader to False to skip use of current importer
1821    def __init__(self,formatName,longFormatName=None,
1822                 extensionlist=[],strictExtension=False,):
1823        self.formatName = formatName # short string naming file type
1824        if longFormatName: # longer string naming file type
1825            self.longFormatName = longFormatName
1826        else:
1827            self.longFormatName = formatName
1828        # define extensions that are allowed for the file type
1829        # for windows, remove any extensions that are duplicate, as case is ignored
1830        if sys.platform == 'windows' and extensionlist:
1831            extensionlist = list(set([s.lower() for s in extensionlist]))
1832        self.extensionlist = extensionlist
1833        # If strictExtension is True, the file will not be read, unless
1834        # the extension matches one in the extensionlist
1835        self.strictExtension = strictExtension
1836        self.errors = ''
1837        self.warnings = ''
1838        self.SciPy = False          #image reader needed scipy
1839        # used for readers that will use multiple passes to read
1840        # more than one data block
1841        self.repeat = False
1842        self.selections = []
1843        self.repeatcount = 0
1844        self.readfilename = '?'
1845        self.scriptable = False
1846        #print 'created',self.__class__
1847
1848    def ReInitialize(self):
1849        'Reinitialize the Reader to initial settings'
1850        self.errors = ''
1851        self.warnings = ''
1852        self.SciPy = False          #image reader needed scipy
1853        self.repeat = False
1854        self.repeatcount = 0
1855        self.readfilename = '?'
1856
1857       
1858#    def Reader(self, filename, filepointer, ParentFrame=None, **unused):
1859#        '''This method must be supplied in the child class to read the file.
1860#        if the read fails either return False or raise an Exception
1861#        preferably of type ImportException.
1862#        '''
1863#        #start reading
1864#        raise ImportException("Error occurred while...")
1865#        self.errors += "Hint for user on why the error occur
1866#        return False # if an error occurs
1867#        return True # if read OK
1868
1869    def ExtensionValidator(self, filename):
1870        '''This methods checks if the file has the correct extension
1871        Return False if this filename will not be supported by this reader
1872        Return True if the extension matches the list supplied by the reader
1873        Return None if the reader allows un-registered extensions
1874        '''
1875        if filename:
1876            ext = ospath.splitext(filename)[1]
1877            if sys.platform == 'windows': ext = ext.lower()
1878            if ext in self.extensionlist: return True
1879            if self.strictExtension: return False
1880        return None
1881
1882    def ContentsValidator(self, filepointer):
1883        '''This routine will attempt to determine if the file can be read
1884        with the current format.
1885        This will typically be overridden with a method that
1886        takes a quick scan of [some of]
1887        the file contents to do a "sanity" check if the file
1888        appears to match the selected format.
1889        '''
1890        #filepointer.seek(0) # rewind the file pointer
1891        return True
1892
1893    def CIFValidator(self, filepointer):
1894        '''A :meth:`ContentsValidator` for use to validate CIF files.
1895        '''
1896        for i,l in enumerate(filepointer):
1897            if i >= 1000: return True
1898            '''Encountered only blank lines or comments in first 1000
1899            lines. This is unlikely, but assume it is CIF anyway, since we are
1900            even less likely to find a file with nothing but hashes and
1901            blank lines'''
1902            line = l.strip()
1903            if len(line) == 0: # ignore blank lines
1904                continue 
1905            elif line.startswith('#'): # ignore comments
1906                continue 
1907            elif line.startswith('data_'): # on the right track, accept this file
1908                return True
1909            else: # found something invalid
1910                self.errors = 'line '+str(i+1)+' contains unexpected data:\n'
1911                if all([ord(c) < 128 and ord(c) != 0 for c in str(l)]): # show only if ASCII
1912                    self.errors += '  '+str(l)
1913                else: 
1914                    self.errors += '  (binary)'
1915                self.errors += '\n  Note: a CIF should only have blank lines or comments before'
1916                self.errors += '\n        a data_ statement begins a block.'
1917                return False 
1918
1919######################################################################
1920class ImportPhase(ImportBaseclass):
1921    '''Defines a base class for the reading of files with coordinates
1922
1923    Objects constructed that subclass this (in import/G2phase_*.py etc.) will be used
1924    in :meth:`GSASIIdataGUI.GSASII.OnImportPhase`.
1925    See :ref:`Writing a Import Routine<Import_Routines>`
1926    for an explanation on how to use this class.
1927
1928    '''
1929    def __init__(self,formatName,longFormatName=None,extensionlist=[],
1930        strictExtension=False,):
1931        # call parent __init__
1932        ImportBaseclass.__init__(self,formatName,longFormatName,
1933            extensionlist,strictExtension)
1934        self.Phase = None # a phase must be created with G2IO.SetNewPhase in the Reader
1935        self.Constraints = None
1936
1937######################################################################
1938class ImportStructFactor(ImportBaseclass):
1939    '''Defines a base class for the reading of files with tables
1940    of structure factors.
1941
1942    Structure factors are read with a call to :meth:`GSASIIdataGUI.GSASII.OnImportSfact`
1943    which in turn calls :meth:`GSASIIdataGUI.GSASII.OnImportGeneric`, which calls
1944    methods :meth:`ExtensionValidator`, :meth:`ContentsValidator` and
1945    :meth:`Reader`.
1946
1947    See :ref:`Writing a Import Routine<Import_Routines>`
1948    for an explanation on how to use import classes in general. The specifics
1949    for reading a structure factor histogram require that
1950    the ``Reader()`` routine in the import
1951    class need to do only a few things: It
1952    should load :attr:`RefDict` item ``'RefList'`` with the reflection list,
1953    and set :attr:`Parameters` with the instrument parameters
1954    (initialized with :meth:`InitParameters` and set with :meth:`UpdateParameters`).
1955    '''
1956    def __init__(self,formatName,longFormatName=None,extensionlist=[],
1957        strictExtension=False,):
1958        ImportBaseclass.__init__(self,formatName,longFormatName,
1959            extensionlist,strictExtension)
1960
1961        # define contents of Structure Factor entry
1962        self.Parameters = []
1963        'self.Parameters is a list with two dicts for data parameter settings'
1964        self.InitParameters()
1965        self.RefDict = {'RefList':[],'FF':{},'Super':0}
1966        self.Banks = []             #for multi bank data (usually TOF)
1967        '''self.RefDict is a dict containing the reflection information, as read from the file.
1968        Item 'RefList' contains the reflection information. See the
1969        :ref:`Single Crystal Reflection Data Structure<XtalRefl_table>`
1970        for the contents of each row. Dict element 'FF'
1971        contains the form factor values for each element type; if this entry
1972        is left as initialized (an empty list) it will be initialized as needed later.
1973        '''
1974    def ReInitialize(self):
1975        'Reinitialize the Reader to initial settings'
1976        ImportBaseclass.ReInitialize(self)
1977        self.InitParameters()
1978        self.Banks = []             #for multi bank data (usually TOF)
1979        self.RefDict = {'RefList':[],'FF':{},'Super':0}
1980       
1981    def InitParameters(self):
1982        'initialize the instrument parameters structure'
1983        Lambda = 0.70926
1984        HistType = 'SXC'
1985        self.Parameters = [{'Type':[HistType,HistType], # create the structure
1986                            'Lam':[Lambda,Lambda]
1987                            }, {}]
1988        'Parameters is a list with two dicts for data parameter settings'
1989
1990    def UpdateParameters(self,Type=None,Wave=None):
1991        'Revise the instrument parameters'
1992        if Type is not None:
1993            self.Parameters[0]['Type'] = [Type,Type]
1994        if Wave is not None:
1995            self.Parameters[0]['Lam'] = [Wave,Wave]           
1996                       
1997######################################################################
1998class ImportPowderData(ImportBaseclass):
1999    '''Defines a base class for the reading of files with powder data.
2000
2001    Objects constructed that subclass this (in import/G2pwd_*.py etc.) will be used
2002    in :meth:`GSASIIdataGUI.GSASII.OnImportPowder`.
2003    See :ref:`Writing a Import Routine<Import_Routines>`
2004    for an explanation on how to use this class.
2005    '''
2006    def __init__(self,formatName,longFormatName=None,
2007        extensionlist=[],strictExtension=False,):
2008        ImportBaseclass.__init__(self,formatName,longFormatName,
2009            extensionlist,strictExtension)
2010        self.clockWd = None  # used in TOF
2011        self.ReInitialize()
2012       
2013    def ReInitialize(self):
2014        'Reinitialize the Reader to initial settings'
2015        ImportBaseclass.ReInitialize(self)
2016        self.powderentry = ['',None,None] #  (filename,Pos,Bank)
2017        self.powderdata = [] # Powder dataset
2018        '''A powder data set is a list with items [x,y,w,yc,yb,yd]:
2019                np.array(x), # x-axis values
2020                np.array(y), # powder pattern intensities
2021                np.array(w), # 1/sig(intensity)^2 values (weights)
2022                np.array(yc), # calc. intensities (zero)
2023                np.array(yb), # calc. background (zero)
2024                np.array(yd), # obs-calc profiles
2025        '''                           
2026        self.comments = []
2027        self.idstring = ''
2028        self.Sample = SetDefaultSample() # default sample parameters
2029        self.Controls = {}  # items to be placed in top-level Controls
2030        self.GSAS = None     # used in TOF
2031        self.repeat_instparm = True # Should a parm file be
2032        #                             used for multiple histograms?
2033        self.instparm = None # name hint from file of instparm to use
2034        self.instfile = '' # full path name to instrument parameter file
2035        self.instbank = '' # inst parm bank number
2036        self.instmsg = ''  # a label that gets printed to show
2037                           # where instrument parameters are from
2038        self.numbanks = 1
2039        self.instdict = {} # place items here that will be transferred to the instrument parameters
2040        self.pwdparms = {} # place parameters that are transferred directly to the tree
2041                           # here (typically from an existing GPX file)
2042######################################################################
2043class ImportSmallAngleData(ImportBaseclass):
2044    '''Defines a base class for the reading of files with small angle data.
2045    See :ref:`Writing a Import Routine<Import_Routines>`
2046    for an explanation on how to use this class.
2047    '''
2048    def __init__(self,formatName,longFormatName=None,extensionlist=[],
2049        strictExtension=False,):
2050           
2051        ImportBaseclass.__init__(self,formatName,longFormatName,extensionlist,
2052            strictExtension)
2053        self.ReInitialize()
2054       
2055    def ReInitialize(self):
2056        'Reinitialize the Reader to initial settings'
2057        ImportBaseclass.ReInitialize(self)
2058        self.smallangleentry = ['',None,None] #  (filename,Pos,Bank)
2059        self.smallangledata = [] # SASD dataset
2060        '''A small angle data set is a list with items [x,y,w,yc,yd]:
2061                np.array(x), # x-axis values
2062                np.array(y), # powder pattern intensities
2063                np.array(w), # 1/sig(intensity)^2 values (weights)
2064                np.array(yc), # calc. intensities (zero)
2065                np.array(yd), # obs-calc profiles
2066                np.array(yb), # preset bkg
2067        '''                           
2068        self.comments = []
2069        self.idstring = ''
2070        self.Sample = SetDefaultSample()
2071        self.GSAS = None     # used in TOF
2072        self.clockWd = None  # used in TOF
2073        self.numbanks = 1
2074        self.instdict = {} # place items here that will be transferred to the instrument parameters
2075
2076######################################################################
2077class ImportReflectometryData(ImportBaseclass):
2078    '''Defines a base class for the reading of files with reflectometry data.
2079    See :ref:`Writing a Import Routine<Import_Routines>`
2080    for an explanation on how to use this class.
2081    '''
2082    def __init__(self,formatName,longFormatName=None,extensionlist=[],
2083        strictExtension=False,):
2084           
2085        ImportBaseclass.__init__(self,formatName,longFormatName,extensionlist,
2086            strictExtension)
2087        self.ReInitialize()
2088       
2089    def ReInitialize(self):
2090        'Reinitialize the Reader to initial settings'
2091        ImportBaseclass.ReInitialize(self)
2092        self.reflectometryentry = ['',None,None] #  (filename,Pos,Bank)
2093        self.reflectometrydata = [] # SASD dataset
2094        '''A small angle data set is a list with items [x,y,w,yc,yd]:
2095                np.array(x), # x-axis values
2096                np.array(y), # powder pattern intensities
2097                np.array(w), # 1/sig(intensity)^2 values (weights)
2098                np.array(yc), # calc. intensities (zero)
2099                np.array(yd), # obs-calc profiles
2100                np.array(yb), # preset bkg
2101        '''                           
2102        self.comments = []
2103        self.idstring = ''
2104        self.Sample = SetDefaultSample()
2105        self.GSAS = None     # used in TOF
2106        self.clockWd = None  # used in TOF
2107        self.numbanks = 1
2108        self.instdict = {} # place items here that will be transferred to the instrument parameters
2109
2110######################################################################
2111class ImportPDFData(ImportBaseclass):
2112    '''Defines a base class for the reading of files with PDF G(R) data.
2113    See :ref:`Writing a Import Routine<Import_Routines>`
2114    for an explanation on how to use this class.
2115    '''
2116    def __init__(self,formatName,longFormatName=None,extensionlist=[],
2117        strictExtension=False,):
2118           
2119        ImportBaseclass.__init__(self,formatName,longFormatName,extensionlist,
2120            strictExtension)
2121        self.ReInitialize()
2122       
2123    def ReInitialize(self):
2124        'Reinitialize the Reader to initial settings'
2125        ImportBaseclass.ReInitialize(self)
2126        self.pdfentry = ['',None,None] #  (filename,Pos,Bank)
2127        self.pdfdata = [] # PDF G(R) dataset
2128        '''A pdf g(r) data set is a list with items [x,y]:
2129                np.array(x), # r-axis values
2130                np.array(y), # pdf g(r)
2131        '''                           
2132        self.comments = []
2133        self.idstring = ''
2134        self.numbanks = 1
2135
2136######################################################################
2137class ImportImage(ImportBaseclass):
2138    '''Defines a base class for the reading of images
2139
2140    Images are read in only these places:
2141   
2142      * Initial reading is typically done from a menu item
2143        with a call to :meth:`GSASIIdataGUI.GSASII.OnImportImage`
2144        which in turn calls :meth:`GSASIIdataGUI.GSASII.OnImportGeneric`. That calls
2145        methods :meth:`ExtensionValidator`, :meth:`ContentsValidator` and
2146        :meth:`Reader`. This returns a list of reader objects for each read image.
2147
2148      * Images are read alternatively in :func:`GSASIIIO.ReadImages`, which puts image info
2149        directly into the data tree.
2150
2151      * Images are reloaded with :func:`GSASIIIO.GetImageData`.
2152
2153    .. _Image_import_routines:
2154
2155    When reading an image, the ``Reader()`` routine in the ImportImage class
2156    should set:
2157   
2158      * :attr:`Comments`: a list of strings (str),
2159      * :attr:`Npix`: the number of pixels in the image (int),
2160      * :attr:`Image`: the actual image as a numpy array (np.array)
2161      * :attr:`Data`: a dict defining image parameters (dict). Within this dict the following
2162        data items are needed:
2163       
2164         * 'pixelSize': size of each pixel in microns (such as ``[200,200]``.
2165         * 'wavelength': wavelength in Angstoms.
2166         * 'distance': distance of detector from sample in cm.
2167         * 'center': uncalibrated center of beam on detector (such as ``[204.8,204.8]``.
2168         * 'size': size of image (such as ``[2048,2048]``).
2169         * 'ImageTag': image number or other keyword used to retrieve image from
2170           a multi-image data file (defaults to ``1`` if not specified).
2171         * 'sumfile': holds sum image file name if a sum was produced from a multi image file
2172
2173    optional data items:
2174   
2175      * :attr:`repeat`: set to True if there are additional images to
2176        read in the file, False otherwise
2177      * :attr:`repeatcount`: set to the number of the image.
2178     
2179    Note that the above is initialized with :meth:`InitParameters`.
2180    (Also see :ref:`Writing a Import Routine<Import_Routines>`
2181    for an explanation on how to use import classes in general.)
2182    '''
2183    def __init__(self,formatName,longFormatName=None,extensionlist=[],
2184        strictExtension=False,):
2185        ImportBaseclass.__init__(self,formatName,longFormatName,
2186            extensionlist,strictExtension)
2187        self.InitParameters()
2188       
2189    def ReInitialize(self):
2190        'Reinitialize the Reader to initial settings -- not used at present'
2191        ImportBaseclass.ReInitialize(self)
2192        self.InitParameters()
2193       
2194    def InitParameters(self):
2195        'initialize the instrument parameters structure'
2196        self.Comments = ['No comments']
2197        self.Data = {}
2198        self.Npix = 0
2199        self.Image = None
2200        self.repeat = False
2201        self.repeatcount = 1
2202        self.sumfile = ''
2203
2204    def LoadImage(self,ParentFrame,imagefile,imagetag=None):
2205        '''Optionally, call this after reading in an image to load it into the tree.
2206        This saves time by preventing a reread of the same information.
2207        '''
2208        if ParentFrame:
2209            ParentFrame.ImageZ = self.Image   # store the image for plotting
2210            ParentFrame.oldImagefile = imagefile # save the name of the last image file read
2211            ParentFrame.oldImageTag = imagetag   # save the tag of the last image file read           
2212
2213#################################################################################################
2214# shortcut routines
2215exp = np.exp
2216sind = sin = s = lambda x: np.sin(x*np.pi/180.)
2217cosd = cos = c = lambda x: np.cos(x*np.pi/180.)
2218tand = tan = t = lambda x: np.tan(x*np.pi/180.)
2219sqrt = sq = lambda x: np.sqrt(x)
2220pi = lambda: np.pi
2221class ExpressionObj(object):
2222    '''Defines an object with a user-defined expression, to be used for
2223    secondary fits or restraints. Object is created null, but is changed
2224    using :meth:`LoadExpression`. This contains only the minimum
2225    information that needs to be stored to save and load the expression
2226    and how it is mapped to GSAS-II variables.
2227    '''
2228    def __init__(self):
2229        self.expression = ''
2230        'The expression as a text string'
2231        self.assgnVars = {}
2232        '''A dict where keys are label names in the expression mapping to a GSAS-II
2233        variable. The value a G2 variable name.
2234        Note that the G2 variable name may contain a wild-card and correspond to
2235        multiple values.
2236        '''
2237        self.freeVars = {}
2238        '''A dict where keys are label names in the expression mapping to a free
2239        parameter. The value is a list with:
2240
2241         * a name assigned to the parameter
2242         * a value for to the parameter and
2243         * a flag to determine if the variable is refined.
2244        ''' 
2245        self.depVar = None
2246
2247        self.lastError = ('','')
2248        '''Shows last encountered error in processing expression
2249        (list of 1-3 str values)'''
2250
2251        self.distance_dict  = None  # to be used for defining atom phase/symmetry info
2252        self.distance_atoms = None  # to be used for defining atom distances
2253
2254    def LoadExpression(self,expr,exprVarLst,varSelect,varName,varValue,varRefflag):
2255        '''Load the expression and associated settings into the object. Raises
2256        an exception if the expression is not parsed, if not all functions
2257        are defined or if not all needed parameter labels in the expression
2258        are defined.
2259
2260        This will not test if the variable referenced in these definitions
2261        are actually in the parameter dictionary. This is checked when the
2262        computation for the expression is done in :meth:`SetupCalc`.
2263       
2264        :param str expr: the expression
2265        :param list exprVarLst: parameter labels found in the expression
2266        :param dict varSelect: this will be 0 for Free parameters
2267          and non-zero for expression labels linked to G2 variables.
2268        :param dict varName: Defines a name (str) associated with each free parameter
2269        :param dict varValue: Defines a value (float) associated with each free parameter
2270        :param dict varRefflag: Defines a refinement flag (bool)
2271          associated with each free parameter
2272        '''
2273        self.expression = expr
2274        self.compiledExpr = None
2275        self.freeVars = {}
2276        self.assgnVars = {}
2277        for v in exprVarLst:
2278            if varSelect[v] == 0:
2279                self.freeVars[v] = [
2280                    varName.get(v),
2281                    varValue.get(v),
2282                    varRefflag.get(v),
2283                    ]
2284            else:
2285                self.assgnVars[v] = varName[v]
2286        self.CheckVars()
2287
2288    def EditExpression(self,exprVarLst,varSelect,varName,varValue,varRefflag):
2289        '''Load the expression and associated settings from the object into
2290        arrays used for editing.
2291
2292        :param list exprVarLst: parameter labels found in the expression
2293        :param dict varSelect: this will be 0 for Free parameters
2294          and non-zero for expression labels linked to G2 variables.
2295        :param dict varName: Defines a name (str) associated with each free parameter
2296        :param dict varValue: Defines a value (float) associated with each free parameter
2297        :param dict varRefflag: Defines a refinement flag (bool)
2298          associated with each free parameter
2299
2300        :returns: the expression as a str
2301        '''
2302        for v in self.freeVars:
2303            varSelect[v] = 0
2304            varName[v] = self.freeVars[v][0]
2305            varValue[v] = self.freeVars[v][1]
2306            varRefflag[v] = self.freeVars[v][2]
2307        for v in self.assgnVars:
2308            varSelect[v] = 1
2309            varName[v] = self.assgnVars[v]
2310        return self.expression
2311
2312    def GetVaried(self):
2313        'Returns the names of the free parameters that will be refined'
2314        return ["::"+self.freeVars[v][0] for v in self.freeVars if self.freeVars[v][2]]
2315
2316    def GetVariedVarVal(self):
2317        'Returns the names and values of the free parameters that will be refined'
2318        return [("::"+self.freeVars[v][0],self.freeVars[v][1]) for v in self.freeVars if self.freeVars[v][2]]
2319
2320    def UpdateVariedVars(self,varyList,values):
2321        'Updates values for the free parameters (after a refinement); only updates refined vars'
2322        for v in self.freeVars:
2323            if not self.freeVars[v][2]: continue
2324            if "::"+self.freeVars[v][0] not in varyList: continue
2325            indx = varyList.index("::"+self.freeVars[v][0])
2326            self.freeVars[v][1] = values[indx]
2327
2328    def GetIndependentVars(self):
2329        'Returns the names of the required independent parameters used in expression'
2330        return [self.assgnVars[v] for v in self.assgnVars]
2331
2332    def CheckVars(self):
2333        '''Check that the expression can be parsed, all functions are
2334        defined and that input loaded into the object is internally
2335        consistent. If not an Exception is raised.
2336
2337        :returns: a dict with references to packages needed to
2338          find functions referenced in the expression.
2339        '''
2340        ret = self.ParseExpression(self.expression)
2341        if not ret:
2342            raise Exception("Expression parse error")
2343        exprLblList,fxnpkgdict = ret
2344        # check each var used in expression is defined
2345        defined = self.assgnVars.keys() + self.freeVars.keys()
2346        notfound = []
2347        for var in exprLblList:
2348            if var not in defined:
2349                notfound.append(var)
2350        if notfound:
2351            msg = 'Not all variables defined'
2352            msg1 = 'The following variables were not defined: '
2353            msg2 = ''
2354            for var in notfound:
2355                if msg: msg += ', '
2356                msg += var
2357            self.lastError = (msg1,'  '+msg2)
2358            raise Exception(msg)
2359        return fxnpkgdict
2360
2361    def ParseExpression(self,expr):
2362        '''Parse an expression and return a dict of called functions and
2363        the variables used in the expression. Returns None in case an error
2364        is encountered. If packages are referenced in functions, they are loaded
2365        and the functions are looked up into the modules global
2366        workspace.
2367       
2368        Note that no changes are made to the object other than
2369        saving an error message, so that this can be used for testing prior
2370        to the save.
2371
2372        :returns: a list of used variables
2373        '''
2374        self.lastError = ('','')
2375        import ast
2376        def FindFunction(f):
2377            '''Find the object corresponding to function f
2378            :param str f: a function name such as 'numpy.exp'
2379            :returns: (pkgdict,pkgobj) where pkgdict contains a dict
2380              that defines the package location(s) and where pkgobj
2381              defines the object associated with the function.
2382              If the function is not found, pkgobj is None.
2383            '''
2384            df = f.split('.')
2385            pkgdict = {}
2386            # no listed package, try in current namespace
2387            if len(df) == 1: 
2388                try:
2389                    fxnobj = eval(f)
2390                    return pkgdict,fxnobj
2391                except (AttributeError, NameError):
2392                    return None,None
2393            else:
2394                try:
2395                    fxnobj = eval(f)
2396                    pkgdict[df[0]] = eval(df[0])
2397                    return pkgdict,fxnobj
2398                except (AttributeError, NameError):
2399                    pass
2400            # includes a package, lets try to load the packages
2401            pkgname = ''
2402            path = sys.path
2403            for pkg in f.split('.')[:-1]: # if needed, descend down the tree
2404                if pkgname:
2405                    pkgname += '.' + pkg
2406                else:
2407                    pkgname = pkg
2408                fp = None
2409                try:
2410                    fp, fppath,desc = imp.find_module(pkg,path)
2411                    pkgobj = imp.load_module(pkg,fp,fppath,desc)
2412                    pkgdict[pkgname] = pkgobj
2413                    path = [fppath]
2414                except Exception as msg:
2415                    print('load of '+pkgname+' failed with error='+str(msg))
2416                    return {},None
2417                finally:
2418                    if fp: fp.close()
2419                try:
2420                    #print 'before',pkgdict.keys()
2421                    fxnobj = eval(f,globals(),pkgdict)
2422                    #print 'after 1',pkgdict.keys()
2423                    #fxnobj = eval(f,pkgdict)
2424                    #print 'after 2',pkgdict.keys()
2425                    return pkgdict,fxnobj
2426                except:
2427                    continue
2428            return None # not found
2429        def ASTtransverse(node,fxn=False):
2430            '''Transverse a AST-parsed expresson, compiling a list of variables
2431            referenced in the expression. This routine is used recursively.
2432
2433            :returns: varlist,fxnlist where
2434              varlist is a list of referenced variable names and
2435              fxnlist is a list of used functions
2436            '''
2437            varlist = []
2438            fxnlist = []
2439            if isinstance(node, list):
2440                for b in node:
2441                    v,f = ASTtransverse(b,fxn)
2442                    varlist += v
2443                    fxnlist += f
2444            elif isinstance(node, ast.AST):
2445                for a, b in ast.iter_fields(node):
2446                    if isinstance(b, ast.AST):
2447                        if a == 'func': 
2448                            fxnlist += ['.'.join(ASTtransverse(b,True)[0])]
2449                            continue
2450                        v,f = ASTtransverse(b,fxn)
2451                        varlist += v
2452                        fxnlist += f
2453                    elif isinstance(b, list):
2454                        v,f = ASTtransverse(b,fxn)
2455                        varlist += v
2456                        fxnlist += f
2457                    elif node.__class__.__name__ == "Name":
2458                        varlist += [b]
2459                    elif fxn and node.__class__.__name__ == "Attribute":
2460                        varlist += [b]
2461            return varlist,fxnlist
2462        try:
2463            exprast = ast.parse(expr)
2464        except SyntaxError:
2465            s = ''
2466            import traceback
2467            for i in traceback.format_exc().splitlines()[-3:-1]:
2468                if s: s += "\n"
2469                s += str(i)
2470            self.lastError = ("Error parsing expression:",s)
2471            return
2472        # find the variables & functions
2473        v,f = ASTtransverse(exprast)
2474        varlist = sorted(list(set(v)))
2475        fxnlist = list(set(f))
2476        pkgdict = {}
2477        # check the functions are defined
2478        for fxn in fxnlist:
2479            fxndict,fxnobj = FindFunction(fxn)
2480            if not fxnobj:
2481                self.lastError = ("Error: Invalid function",fxn,
2482                                  "is not defined")
2483                return
2484            if not hasattr(fxnobj,'__call__'):
2485                self.lastError = ("Error: Not a function.",fxn,
2486                                  "cannot be called as a function")
2487                return
2488            pkgdict.update(fxndict)
2489        return varlist,pkgdict
2490
2491    def GetDepVar(self):
2492        'return the dependent variable, or None'
2493        return self.depVar
2494
2495    def SetDepVar(self,var):
2496        'Set the dependent variable, if used'
2497        self.depVar = var
2498#==========================================================================
2499class ExpressionCalcObj(object):
2500    '''An object used to evaluate an expression from a :class:`ExpressionObj`
2501    object.
2502   
2503    :param ExpressionObj exprObj: a :class:`~ExpressionObj` expression object with
2504      an expression string and mappings for the parameter labels in that object.
2505    '''
2506    def __init__(self,exprObj):
2507        self.eObj = exprObj
2508        'The expression and mappings; a :class:`ExpressionObj` object'
2509        self.compiledExpr = None
2510        'The expression as compiled byte-code'
2511        self.exprDict = {}
2512        '''dict that defines values for labels used in expression and packages
2513        referenced by functions
2514        '''
2515        self.lblLookup = {}
2516        '''Lookup table that specifies the expression label name that is
2517        tied to a particular GSAS-II parameters in the parmDict.
2518        '''
2519        self.fxnpkgdict = {}
2520        '''a dict with references to packages needed to
2521        find functions referenced in the expression.
2522        '''
2523        self.varLookup = {}
2524        '''Lookup table that specifies the GSAS-II variable(s)
2525        indexed by the expression label name. (Used for only for diagnostics
2526        not evaluation of expression.)
2527        '''
2528        self.su = None
2529        '''Standard error evaluation where supplied by the evaluator
2530        '''
2531        # Patch: for old-style expressions with a (now removed step size)
2532        for v in self.eObj.assgnVars:
2533            if not isinstance(self.eObj.assgnVars[v], basestring):
2534                self.eObj.assgnVars[v] = self.eObj.assgnVars[v][0]
2535        self.parmDict = {}
2536        '''A copy of the parameter dictionary, for distance and angle computation
2537        '''
2538
2539    def SetupCalc(self,parmDict):
2540        '''Do all preparations to use the expression for computation.
2541        Adds the free parameter values to the parameter dict (parmDict).
2542        '''
2543        if self.eObj.expression.startswith('Dist') or self.eObj.expression.startswith('Angle'):
2544            return
2545        self.fxnpkgdict = self.eObj.CheckVars()
2546        # all is OK, compile the expression
2547        self.compiledExpr = compile(self.eObj.expression,'','eval')
2548
2549        # look at first value in parmDict to determine its type
2550        parmsInList = True
2551        for key in parmDict:
2552            val = parmDict[key]
2553            if isinstance(val, basestring):
2554                parmsInList = False
2555                break
2556            try: # check if values are in lists
2557                val = parmDict[key][0]
2558            except (TypeError,IndexError):
2559                parmsInList = False
2560            break
2561           
2562        # set up the dicts needed to speed computations
2563        self.exprDict = {}
2564        self.lblLookup = {}
2565        self.varLookup = {}
2566        for v in self.eObj.freeVars:
2567            varname = self.eObj.freeVars[v][0]
2568            varname = "::" + varname.lstrip(':').replace(' ','_').replace(':',';')
2569            self.lblLookup[varname] = v
2570            self.varLookup[v] = varname
2571            if parmsInList:
2572                parmDict[varname] = [self.eObj.freeVars[v][1],self.eObj.freeVars[v][2]]
2573            else:
2574                parmDict[varname] = self.eObj.freeVars[v][1]
2575            self.exprDict[v] = self.eObj.freeVars[v][1]
2576        for v in self.eObj.assgnVars:
2577            varname = self.eObj.assgnVars[v]
2578            if '*' in varname:
2579                varlist = LookupWildCard(varname,parmDict.keys())
2580                if len(varlist) == 0:
2581                    raise Exception,"No variables match "+str(v)
2582                for var in varlist:
2583                    self.lblLookup[var] = v
2584                if parmsInList:
2585                    self.exprDict[v] = np.array([parmDict[var][0] for var in varlist])
2586                else:
2587                    self.exprDict[v] = np.array([parmDict[var] for var in varlist])
2588                self.varLookup[v] = [var for var in varlist]
2589            elif varname in parmDict:
2590                self.lblLookup[varname] = v
2591                self.varLookup[v] = varname
2592                if parmsInList:
2593                    self.exprDict[v] = parmDict[varname][0]
2594                else:
2595                    self.exprDict[v] = parmDict[varname]
2596            else:
2597                self.exprDict[v] = None
2598#                raise Exception,"No value for variable "+str(v)
2599        self.exprDict.update(self.fxnpkgdict)
2600
2601    def UpdateVars(self,varList,valList):
2602        '''Update the dict for the expression with a set of values
2603        :param list varList: a list of variable names
2604        :param list valList: a list of corresponding values
2605        '''
2606        for var,val in zip(varList,valList):
2607            self.exprDict[self.lblLookup.get(var,'undefined: '+var)] = val
2608
2609    def UpdateDict(self,parmDict):
2610        '''Update the dict for the expression with values in a dict
2611        :param list parmDict: a dict of values some of which may be in use here
2612        '''
2613        if self.eObj.expression.startswith('Dist') or self.eObj.expression.startswith('Angle'):
2614            self.parmDict = parmDict
2615            return
2616        for var in parmDict:
2617            if var in self.lblLookup:
2618                self.exprDict[self.lblLookup[var]] = parmDict[var]
2619           
2620    def EvalExpression(self):
2621        '''Evaluate an expression. Note that the expression
2622        and mapping are taken from the :class:`ExpressionObj` expression object
2623        and the parameter values were specified in :meth:`SetupCalc`.
2624        :returns: a single value for the expression. If parameter
2625        values are arrays (for example, from wild-carded variable names),
2626        the sum of the resulting expression is returned.
2627
2628        For example, if the expression is ``'A*B'``,
2629        where A is 2.0 and B maps to ``'1::Afrac:*'``, which evaluates to::
2630
2631        [0.5, 1, 0.5]
2632
2633        then the result will be ``4.0``.
2634        '''
2635        self.su = None
2636        if self.eObj.expression.startswith('Dist'):
2637#            GSASIIpath.IPyBreak()
2638            dist = G2mth.CalcDist(self.eObj.distance_dict, self.eObj.distance_atoms, self.parmDict)
2639            return dist
2640        elif self.eObj.expression.startswith('Angle'):
2641            angle = G2mth.CalcAngle(self.eObj.angle_dict, self.eObj.angle_atoms, self.parmDict)
2642            return angle
2643        if self.compiledExpr is None:
2644            raise Exception,"EvalExpression called before SetupCalc"
2645        try:
2646            val = eval(self.compiledExpr,globals(),self.exprDict)
2647        except TypeError:
2648            val = None
2649        if not np.isscalar(val):
2650            val = np.sum(val)
2651        return val
2652       
2653class G2Exception(Exception):
2654    def __init__(self,msg):
2655        self.msg = msg
2656    def __str__(self):
2657        return repr(self.msg)
2658
2659def HowDidIgetHere(wherecalledonly=False):
2660    '''Show a traceback with calls that brought us to the current location.
2661    Used for debugging.
2662    '''
2663    import traceback
2664    if wherecalledonly:
2665        i = traceback.format_list(traceback.extract_stack()[:-1])[-2]
2666        print(i.strip().rstrip())
2667    else:
2668        print 70*'*'   
2669        for i in traceback.format_list(traceback.extract_stack()[:-1]): print(i.strip().rstrip())
2670        print 70*'*'   
2671               
2672def CreatePDFitems(G2frame,PWDRtree,ElList,Qlimits,numAtm=1,FltBkg=0,PDFnames=[]):
2673    '''Create and initialize a new set of PDF tree entries
2674
2675    :param Frame G2frame: main GSAS-II tree frame object
2676    :param str PWDRtree: name of PWDR to be used to create PDF item
2677    :param dict ElList: data structure with composition
2678    :param list Qlimits: Q limits to be used for computing the PDF
2679    :param float numAtm: no. atom in chemical formula
2680    :param float FltBkg: flat background value
2681    :param list PDFnames: previously used PDF names
2682   
2683    :returns: the Id of the newly created PDF entry
2684    '''
2685    PDFname = 'PDF '+PWDRtree[4:] # this places two spaces after PDF, which is needed is some places
2686    if PDFname in PDFnames:
2687        print('Skipping, entry already exists: '+PDFname)
2688        return None
2689    #PDFname = MakeUniqueLabel(PDFname,PDFnames)
2690    Id = G2frame.GPXtree.AppendItem(parent=G2frame.root,text=PDFname)
2691    Data = {
2692        'Sample':{'Name':PWDRtree,'Mult':1.0},
2693        'Sample Bkg.':{'Name':'','Mult':-1.0,'Refine':False},
2694        'Container':{'Name':'','Mult':-1.0,'Refine':False},
2695        'Container Bkg.':{'Name':'','Mult':-1.0},'ElList':ElList,
2696        'Geometry':'Cylinder','Diam':1.0,'Pack':0.50,'Form Vol':10.0*numAtm,'Flat Bkg':FltBkg,
2697        'DetType':'Area detector','ObliqCoeff':0.2,'Ruland':0.025,'QScaleLim':Qlimits,
2698        'Lorch':False,'BackRatio':0.0,'Rmax':100.,'noRing':False,'IofQmin':1.0,
2699        'I(Q)':[],'S(Q)':[],'F(Q)':[],'G(R)':[]}
2700    G2frame.GPXtree.SetItemPyData(G2frame.GPXtree.AppendItem(Id,text='PDF Controls'),Data)
2701    G2frame.GPXtree.SetItemPyData(G2frame.GPXtree.AppendItem(Id,text='PDF Peaks'),
2702        {'Limits':[1.,5.],'Background':[2,[0.,-0.2*np.pi],False],'Peaks':[]})
2703    return Id       
2704
2705
2706if __name__ == "__main__":
2707    # test equation evaluation
2708    def showEQ(calcobj):
2709        print 50*'='
2710        print calcobj.eObj.expression,'=',calcobj.EvalExpression()
2711        for v in sorted(calcobj.varLookup):
2712            print "  ",v,'=',calcobj.exprDict[v],'=',calcobj.varLookup[v]
2713        # print '  Derivatives'
2714        # for v in calcobj.derivStep.keys():
2715        #     print '    d(Expr)/d('+v+') =',calcobj.EvalDeriv(v)
2716
2717    obj = ExpressionObj()
2718
2719    obj.expression = "A*np.exp(B)"
2720    obj.assgnVars =  {'B': '0::Afrac:1'}
2721    obj.freeVars =  {'A': [u'A', 0.5, True]}
2722    #obj.CheckVars()
2723    calcobj = ExpressionCalcObj(obj)
2724
2725    obj1 = ExpressionObj()
2726    obj1.expression = "A*np.exp(B)"
2727    obj1.assgnVars =  {'B': '0::Afrac:*'}
2728    obj1.freeVars =  {'A': [u'Free Prm A', 0.5, True]}
2729    #obj.CheckVars()
2730    calcobj1 = ExpressionCalcObj(obj1)
2731
2732    obj2 = ExpressionObj()
2733    obj2.distance_stuff = np.array([[0,1],[1,-1]])
2734    obj2.expression = "Dist(1,2)"
2735    GSASIIpath.InvokeDebugOpts()
2736    parmDict2 = {'0::Afrac:0':[0.0,True], '0::Afrac:1': [1.0,False]}
2737    calcobj2 = ExpressionCalcObj(obj2)
2738    calcobj2.SetupCalc(parmDict2)
2739    showEQ(calcobj2)
2740   
2741    parmDict1 = {'0::Afrac:0':1.0, '0::Afrac:1': 1.0}
2742    print '\nDict = ',parmDict1
2743    calcobj.SetupCalc(parmDict1)
2744    showEQ(calcobj)
2745    calcobj1.SetupCalc(parmDict1)
2746    showEQ(calcobj1)
2747
2748    parmDict2 = {'0::Afrac:0':[0.0,True], '0::Afrac:1': [1.0,False]}
2749    print 'Dict = ',parmDict2
2750    calcobj.SetupCalc(parmDict2)
2751    showEQ(calcobj)
2752    calcobj1.SetupCalc(parmDict2)
2753    showEQ(calcobj1)
2754    calcobj2.SetupCalc(parmDict2)
2755    showEQ(calcobj2)
Note: See TracBrowser for help on using the repository browser.