source: trunk/GSASIIobj.py @ 1057

Last change on this file since 1057 was 1057, checked in by toby, 8 years ago

reorg UserCalibrants? into ImageCalibrants?; Add trim (whitespace) routine; fix InstName? bug in GetHistogramPhaseData?; rework ValidatedTxtCtrl? for CIF; ScrolledMultiEditor? now resets after Cancel; CIF export progress

  • Property svn:eol-style set to native
File size: 35.6 KB
Line 
1# TODO: change this to assemble the look-up tables of atoms, phases and hists from the tree
2# and then save/unsave those values in __init__ & __str__, etc.
3
4'''
5*GSASIIobj: Data objects*
6=========================
7
8This module defines and/or documents the data structures used in GSAS-II.
9
10
11Constraints Tree Item
12----------------------
13
14.. _Constraints_table:
15
16.. index::
17   single: Constraints object description
18   single: Data object descriptions; Constraints
19
20Constraints are stored in a dict, separated into groups.
21Note that parameter are named in the following pattern,
22p:h:<var>:n, where p is the phase number, h is the histogram number
23<var> is a variable name and n is the parameter number.
24If a parameter does not depend on a histogram or phase or is unnumbered, that
25number is omitted.
26Note that the contents of each dict item is a List where each element in the
27list is a :ref:`constraint definition objects <Constraint_definitions_table>`.
28
29The keys in the Constraints dict are:
30
31.. tabularcolumns:: |l|p{4.5in}|
32
33==========  ====================================================
34  key         explanation
35==========  ====================================================
36Hist        This specifies a list of constraints on
37            histogram-related parameters,
38            which will be of form :h:<var>:n.
39HAP         This specifies a list of constraints on parameters
40            that are defined for every histogram in each phase
41            and are of form p:h:<var>:n.           
42Phase       This specifies a list of constraints on phase
43            parameters,
44            which will be of form p::<var>:n.
45Global      This specifies a list of constraints on parameters
46            that are not tied to a histogram or phase and
47            are of form ::<var>:n
48==========  ====================================================
49
50.. _Constraint_definitions_table:
51
52.. index::
53   single: Constraint definition object description
54   single: Data object descriptions; Constraint Definition
55
56Each constraint is defined as a list using a series of terms of form
57
58::
59
60[[<mult1>, <var1>], [<mult2>, <var2>],..., <fixed val>, <vary flag>, <cons type>]
61
62Where the variable pair list item containing two values [<mult>, <var>],
63
64  * <mult> is a multiplier for the constraint (float)
65  * <var> is the name of the variable (str) (or to be implemented a :class:`VarName` object.)
66
67Note that the last three items in the list play a special role:
68
69 * <fixed val> is the fixed value for a constraint equation or is None
70 * <vary flag> is True, False or None and is intended to use to indicate if new variables
71   should be refined.
72 * <cons type> is one of four letters, 'e', 'c', 'h', 'f' that determines the type of constraint.
73
74    * 'e' defines a set of equivalent variables. Only the first variable is refined (if the
75      appropriate refine flag is set) and and all other equivalent variables in the list
76      are generated from that variable. The vary flag for those variables is ignored.
77    * 'c' defines a constraint equation of form, :math:`m_1 \\times var_1 + m_2 \\times var_2 + ... = c`
78    * 'h' defines a variable to hold (not vary). Any variable on this list is not varied, even if its refinement
79      flag is set. This is of particular value when needing to hold one or more variables in a set such as
80      the reciprocal metric tensor or anisotropic displacement parameter.
81    * 'f' defines a relationship to define a new variable according to relationship
82      :math:`newvar = m_1 \\times var_1 + m_2 \\times var_2 + ...`
83
84Covariance Tree Item
85--------------------
86
87.. _Covariance_table:
88
89.. index::
90   single: Covariance description
91   single: Data object descriptions; Covariance
92
93The Covariance tree item has results from the last least-squares run. They
94are stored in a dict with these keys:
95
96.. tabularcolumns:: |l|l|p{4in}|
97
98=============  ===============  ====================================================
99  key            sub-key        explanation
100=============  ===============  ====================================================
101newCellDict    \                dict with lattice parameters computed by
102                                :func:`GSASIIstrMath.GetNewCellParms` (dict)
103title          \                Name of gpx file(?) (str)
104variables      \                Values for all N refined variables
105                                (list of float values, length N,
106                                ordered to match varyList)
107sig            \                Uncertainty values for all N refined variables
108                                (list of float values, length N,
109                                ordered to match varyList)
110varyList       \                List of directly refined variables
111                                (list of str values, length N)
112newAtomDict    \                dict with atom position values computed in
113                                :func:`GSASIIstrMath.ApplyXYZshifts` (dict)
114Rvals          \                R-factors, GOF, Marquardt value for last
115                                refinement cycle (dict)
116\              Nobs             Number of observed data points (int)
117\              Rwp              overall weighted profile R-factor (%, float)
118\              chisq            sum[w*(Iobs-Icalc)**2] for all data
119                                note this is not the reduced chi squared (float)
120\              lamMax           Marquardt value applied to Hessian diagonal
121                                (float)
122\              GOF              The goodness-of-fit, aka square root of
123                                the reduced chi squared. (float)
124covMatrix      \                The (NxN) covVariance matrix (np.array)
125=============  ===============  ====================================================
126
127Phase Tree Items
128----------------
129
130.. _Phase_table:
131
132.. index::
133   single: Phase object description
134   single: Data object descriptions; Phase
135
136Phase information is stored in the GSAS-II data tree as children of the
137Phases item in a dict with keys:
138
139.. tabularcolumns:: |l|l|p{4in}|
140
141==========  ===============  ====================================================
142  key         sub-key        explanation
143==========  ===============  ====================================================
144General         \            Overall information for the phase (dict)
145  \         AtomPtrs         list of four locations to use to pull info
146                             from the atom records (list)
147  \         F000X            x-ray F(000) intensity (float)
148  \         F000N            neutron F(000) intensity (float)
149  \         Mydir            directory of current .gpx file (str)
150  \         MCSA controls    ?
151  \         Cell             List with 7 items: cell refinement flag (bool)
152                             a, b, c, (Angstrom, float)
153                             alpha, beta & gamma (degrees, float)
154  \         Type             for now 'nuclear' (str)
155  \         Map              dict of map parameters
156  \         SH Texture       dict of spherical harmonic preferred orientation
157                             parameters
158  \         Isotope          dict of isotopes for each atom type
159  \         Isotopes         dict of scattering lengths for each isotope
160                             combination for each element in phase 
161  \         Name             phase name (str)
162  \         SGData           Space group details as a :ref:`space group (SGData) object <SGData_table>`
163                             as defined in :func:`GSASIIspc.SpcGroup`.
164  \         Pawley neg wt    Restraint value for negative Pawley intensities
165                             (float)
166  \         Flip             Charge flip controls dict?
167  \         Data plot type   ?
168  \         Mass             Mass of unit cell contents in g/mol
169  \         POhkl            March-Dollase preferred orientation direction
170  \         Z                ?
171  \         vdWRadii         ?
172  \         Color            Colors for atoms (list of (r,b,g) triplets)
173  \         AtomTypes        List of atom types
174  \         AtomMass         List of masses for atoms
175  \         doPawley         Flag for Pawley intensity extraction (bool)
176  \         NoAtoms          Number of atoms per unit cell of each type (dict)
177  \         Pawley dmin      maximum Q (as d-space) to use for Pawley
178                             extraction (float)
179  \         BondRadii        Radius for each atom used to compute
180                             interatomic distances (list of floats)
181  \         AngleRadii       Radius for each atom used to compute
182                             interatomic angles (list of floats)
183ranId           \            unique random number Id for phase (int)
184pId             \            Phase Id number for current project (int).
185Atoms           \            Atoms in phase as a list of lists. The outer list
186                             is for each atom, the inner list contains varying
187                             items depending on the type of phase, see
188                             the :ref:`Atom Records <Atoms_table>` description.
189                             (list of lists)
190Drawing         \            Display parameters (dict)
191\           ballScale        Size of spheres in ball-and-stick display (float)
192\           bondList         dict with bonds
193\           contourLevel     ? (float)
194\           showABC          Flag to show view point triplet (bool). True=show.
195\           viewDir          cartesian viewing direction (np.array with three
196                             elements)
197\           Zclip            clipping distance in A (float)
198\           backColor        background for plot as and R,G,B triplet
199                             (default = [0, 0, 0], black).
200                             (list with three atoms)
201\           selectedAtoms    List of selected atoms (list of int values)
202\           showRigidBodies  Flag to highlight rigid body placement
203\           sizeH            Size ratio for H atoms (float)
204\           bondRadius       Size of binds in A (float)
205\           atomPtrs         ? (list)
206\           viewPoint        list of lists. First item in list is [x,y,z]
207                             in fractional coordinates for the center of
208                             the plot. Second item ?.
209\           showHydrogen     Flag to control plotting of H atoms.
210\           unitCellBox      Flag to control display of the unit cell.
211\           ellipseProb      Probability limit for display of thermal
212                             ellipsoids in % (float).
213\           vdwScale         Multiplier of van der Waals radius for
214                             display of vdW spheres.
215\           Atoms            A list of lists with an entry for each atom
216                             that is plotted.
217\           Zstep            Step to de/increase Z-clip (float)
218\           Quaternion       Viewing quaternion (4 element np.array)
219\           radiusFactor     Distance ratio for searching for bonds. ? Bonds
220                             are located that are within r(Ra+Rb) and (Ra+Rb)/r
221                             where Ra and Rb are the atomic radii.
222\           oldxy            ? (list with two floats)
223\           cameraPos        Viewing position in A for plot (float)
224\           depthFog         ? (bool)
225RBModels        \            Rigid body assignments (note Rigid body definitions
226                             are stored in their own main top-level tree entry.)
227Pawley ref      \            Pawley reflections
228Histograms      \            A dict of dicts. The key for the outer dict is
229                             the histograms tied to this phase. The inner
230                             dict contains the combined phase/histogram
231                             parameters for items such as scale factors,
232                             size and strain parameters. (dict)
233MCSA            \            Monte-Carlo simulated annealing parameters
234==========  ===============  ====================================================
235
236Space Group Objects
237-------------------
238
239.. _SGData_table:
240
241.. index::
242   single: SGData description
243   single: Data object descriptions; SGData
244
245Space groups are interpreted by :func:`GSASIIspc.SpcGroup`
246and the information is placed in a SGdata object,
247which is a dict with these keys:
248
249.. tabularcolumns:: |l|p{4.5in}|
250
251==========  ====================================================
252  key         explanation
253==========  ====================================================
254SpGrp       space group symbol (str)
255Laue        one of the following 14 Laue classes:
256            -1, 2/m, mmm, 4/m, 4/mmm, 3R,
257            3mR, 3, 3m1, 31m, 6/m, 6/mmm, m3, m3m (str)
258SGInv       True if centrosymmetric, False if not (bool)
259SGLatt      Lattice centering type. Will be one of
260            P, A, B, C, I, F, R (str)
261SGUniq      unique axis if monoclinic. Will be
262            a, b, or c for monoclinic space groups.
263            Will be blank for non-monoclinic. (str)
264SGCen       Symmetry cell centering vectors. A (n,3) np.array
265            of centers. Will always have at least one row:
266            ``np.array([[0, 0, 0]])``
267SGOps       symmetry operations as a list of form
268            ``[[M1,T1], [M2,T2],...]``
269            where :math:`M_n` is a 3x3 np.array
270            and :math:`T_n` is a length 3 np.array.
271            Atom coordinates are transformed where the
272            Asymmetric unit coordinates [X is (x,y,z)]
273            are transformed using
274            :math:`X^\prime = M_n*X+T_n`
275SGSys       symmetry unit cell: type one of
276            'triclinic', 'monoclinic', 'orthorhombic',
277            'tetragonal', 'rhombohedral', 'trigonal',
278            'hexagonal', 'cubic' (str)
279SGPolax     Axes for space group polarity. Will be one of
280            '', 'x', 'y', 'x y', 'z', 'x z', 'y z',
281            'xyz'. In the case where axes are arbitrary
282            '111' is used (P 1, and ?).
283==========  ====================================================
284
285Atom Records
286------------
287
288.. _Atoms_table:
289
290.. index::
291   single: Atoms record description
292   single: Data object descriptions; Atoms record
293
294
295If ``phasedict`` points to the phase information in the data tree, then
296atoms are contained in a list of atom records (list) in
297``phasedict['Atoms']``. Also needed to read atom information
298are four pointers, ``cx,ct,cs,cia = phasedict['General']['AtomPtrs']``,
299which define locations in the atom record, as shown below.
300
301.. tabularcolumns:: |l|p{4.5in}|
302
303==============   ====================================================
304location         explanation
305==============   ====================================================
306cx,cx+1,cx+2      the x,y and z coordinates
307cx+3              fractional occupancy (also cs-1)
308ct-1              atom label
309ct                atom type
310ct+1              refinement flags
311cs                site symmetry string
312cs+1              site multiplicity
313cia               ADP flag: Isotropic ('I') or Anisotropic ('A')
314cia+1             Uiso
315cia+2...cia+6     U11, U22, U33, U12, U13, U23
316==============   ====================================================
317
318Powder Diffraction Tree Items
319-----------------------------
320
321.. _Powder_table:
322
323.. index::
324   single: Powder data object description
325   single: Data object descriptions; Powder Data
326
327Every powder diffraction histogram is stored in the GSAS-II data tree
328with a top-level entry named beginning with the string "PWDR ". The
329diffraction data for that information are directly associated with
330that tree item and there are a series of children to that item. The
331routine :func:`~GSASII.GSASII.GetUsedHistogramsAndPhasesfromTree` will
332load this information into a dictionary where the child tree name is
333used as a key, and the information in the main entry is assigned
334a key of ``Data``, as outlined below.
335
336.. tabularcolumns:: |l|l|p{4in}|
337
338======================  ===============  ====================================================
339  key                      sub-key        explanation
340======================  ===============  ====================================================
341Limits                       \            A list of two two element lists, as [[Ld,Hd],[L,H]]
342                                          where L and Ld are the current and default lowest
343                                          two-theta value to be used and
344                                          where H and Hd are the current and default highest
345                                          two-theta value to be used.
346Reflection Lists              \           A dict with an entry for each phase in the
347                                          histogram. The contents of each dict item
348                                          is a list of reflections as described in the
349                                          :ref:`Powder Reflections <PowderRefl_table>`
350                                          description.
351Instrument Parameters         \           A list containing two dicts where the possible
352                                          keys in each dict are listed below. The value
353                                          for each item is a list containing three values:
354                                          the initial value, the current value and a
355                                          refinement flag which can have a value of
356                                          True, False or 0 where 0 indicates a value that
357                                          cannot be refined. The first and second
358                                          values are floats unless otherwise noted.
359                                          Items in the first dict are noted as [1]
360\                         Lam             Specifies a wavelength in Angstroms [1]
361\                         Lam1            Specifies the primary wavelength in
362                                          Angstrom, when an alpha1, alpha2
363                                          source is used [1]
364\                         Lam2            Specifies the secondary wavelength in
365                                          Angstrom, when an alpha1, alpha2
366                                          source is used [1]
367                          I(L2)/I(L1)     Ratio of Lam2 to Lam1 [1]           
368\                         Type            Histogram type (str) [1]:
369                                           * 'PXC' for constant wavelength x-ray
370                                           * 'PNC' for constant wavelength neutron
371                                           * 'PNT' for time of flight neutron
372\                         Zero            Two-theta zero correction in *degrees* [1]
373\                         Azimuth         Azimuthal setting angle for data recorded
374                                          with differing setting angles [1]
375\                         U, V, W         Cagliotti profile coefficients
376                                          for Gaussian instrumental broadening, where the
377                                          FWHM goes as
378                                          :math:`U \\tan^2\\theta + V \\tan\\theta + W` [1]
379\                         X, Y            Cauchy (Lorentzian) instrumental broadening
380                                          coefficients [1]
381\                         SH/L            Variant of the Finger-Cox-Jephcoat asymmetric
382                                          peak broadening ratio. Note that this is the
383                                          average between S/L and H/L where S is
384                                          sample height, H is the slit height and
385                                          L is the goniometer diameter. [1]
386\                         Polariz.        Polarization coefficient. [1]
387wtFactor                      \           A weighting factor to increase or decrease
388                                          the leverage of data in the histogram (float).
389                                          A value of 1.0 weights the data with their
390                                          standard uncertainties and a larger value
391                                          increases the weighting of the data (equivalent
392                                          to decreasing the uncertainties).
393Sample Parameters             \           Specifies a dict with parameters that describe how
394                                          the data were collected, as listed
395                                          below. Refinable parameters are a list containing
396                                          a float and a bool, where the second value
397                                          specifies if the value is refined, otherwise
398                                          the value is a float unless otherwise noted.
399\                         Scale           The histogram scale factor (refinable)
400\                         Absorption      The sample absorption coefficient as
401                                          :math:`\\mu r` where r is the radius
402                                          (refinable).
403\                         DisplaceX,      Sample displacement from goniometer center
404                          DisplaceY       where Y is along the beam direction and
405                                          X is perpendicular. Units are :math:`\\mu m`
406                                          (refinable).
407\                         Phi, Chi,       Goniometer sample setting angles, in degrees.
408                          Omega
409\                         Gonio. radius   Radius of the diffractometer in mm
410\                         InstrName       A name for the instrument, used in preparing
411                                          a CIF (str).
412\                         Force,          Variables that describe how the measurement
413                          Temperature,    was performed. Not used directly in
414                          Humidity,       any computations.
415                          Pressure,
416                          Voltage
417\                         ranId           The random-number Id for the histogram
418                                          (same value as where top-level key is ranId)
419\                         Type            Type of diffraction data, may be 'Debye-Scherrer'
420                                          or 'Bragg-Brentano' (str).
421\                         Diffuse         not in use?
422hId                           \           The number assigned to the histogram when
423                                          the project is loaded or edited (can change)
424ranId                         \           A random number id for the histogram
425                                          that does not change
426Background                    \           The background is stored as a list with where
427                                          the first item in the list is list and the second
428                                          item is a dict. The list contains the background
429                                          function and its coefficients; the dict contains
430                                          Debye diffuse terms and background peaks.
431                                          (TODO: this needs to be expanded.)
432Data                          \           The data consist of a list of 6 np.arrays
433                                          containing in order:
434
435                                           1. the x-postions (two-theta in degrees),
436                                           2. the intensity values (Yobs),
437                                           3. the weights for each Yobs value
438                                           4. the computed intensity values (Ycalc)
439                                           5. the background values
440                                           6. Yobs-Ycalc
441======================  ===============  ====================================================
442
443Powder Reflection Data Structure
444--------------------------------
445
446.. _PowderRefl_table:
447
448.. index::
449   single: Powder reflection object description
450   single: Data object descriptions; Powder Reflections
451   
452For every phase in a histogram, the ``Reflection Lists`` value is a list of
453reflections. The items in that list are documented below.
454
455==========  ====================================================
456  index         explanation
457==========  ====================================================
458 0,1,2       h,k,l (float)
459 3           multiplicity
460 4           d-space, Angstrom
461 5           pos, two-theta
462 6           sig, Gaussian width
463 7           gam, Lorenzian width
464 8           :math:`F_{obs}^2`
465 9           :math:`F_{calc}^2`
466 10          reflection phase, in degrees
467 11          the equivalent reflections as a (m x 3)
468             np.array, where m is 0.5 * multiplicity. Note
469             that Freidel pairs, (-h,-k-,l), are not
470             included.
471 12          phase shift for each of the equivalent
472             reflections as a length (m) array
473 13          intensity correction for reflection, this times
474             :math:`F_{obs}^2` or :math:`F_{calc}^2` gives Iobs or Icalc
475 14          dict with the form factor (f or b) by atom type
476             symbol at the reflection position.
477==========  ====================================================
478
479Single Crystal Tree Items
480-------------------------
481
482.. _Xtal_table:
483
484.. index::
485   single: Single Crytsal data object description
486   single: Data object descriptions; Single crystal data
487
488Every single crystal diffraction histogram is stored in the GSAS-II data tree
489with a top-level entry named beginning with the string "HKLF ". The
490diffraction data for that information are directly associated with
491that tree item and there are a series of children to that item. The
492routine :func:`~GSASII.GSASII.GetUsedHistogramsAndPhasesfromTree` will
493load this information into a dictionary where the child tree name is
494used as a key, and the information in the main entry is assigned
495a key of ``Data``, as outlined below.
496
497.. tabularcolumns:: |l|l|p{4in}|
498
499======================  ===============  ====================================================
500  key                      sub-key        explanation
501======================  ===============  ====================================================
502Data
503                                          A list of lists, where each inner item
504                                          is an individual reflection
505                                          as described in the
506                                          :ref:`Single Crystal Reflections
507                                          <XtalRefl_table>`
508                                          description.
509
510Instrument Parameters         \           A list containing two dicts where the possible
511                                          keys in each dict are listed below. The value
512                                          for most items is a list containing two values:
513                                          the initial value, the current value.
514                                          The first and second
515                                          values are floats unless otherwise noted.
516\                         Lam             Specifies a wavelength in Angstroms (two floats)
517\                         Type            Histogram type (two str values):
518                                           * 'SXC' for constant wavelength x-ray
519                                           * 'SNC' for constant wavelength neutron
520                                           * 'SNT' for time of flight neutron
521\                         InstrName       A name for the instrument, used in preparing
522                                          a CIF (str).
523
524wtFactor                      \           A weighting factor to increase or decrease
525                                          the leverage of data in the histogram (float).
526                                          A value of 1.0 weights the data with their
527                                          standard uncertainties and a larger value
528                                          increases the weighting of the data (equivalent
529                                          to decreasing the uncertainties).
530
531hId                           \           The number assigned to the histogram when
532                                          the project is loaded or edited (can change)
533======================  ===============  ====================================================
534
535Single Crystal Reflection Data Structure
536----------------------------------------
537
538.. _XtalRefl_table:
539
540.. index::
541   single: Single Crystal reflection object description
542   single: Data object descriptions; Single Crystal Reflections
543   
544For every phase in a histogram, the ``Reflection Lists`` value is a list of
545reflections. The items in that list are documented below.
546
547==========  ====================================================
548  index         explanation
549==========  ====================================================
550 0,1,2       h,k,l (float)
551 3           multiplicity
552 4           d-space, Angstrom
553 5           :math:`F_{obs}^2`
554 6           :math:`\sigma(F_{obs}^2)`
555 7           :math:`F_{calc}^2`
556 8           :math:`F_{obs}^2T`
557 9           :math:`F_{calc}^2T`
558 10          reflection phase, in degrees
559 11          the equivalent reflections as a (m x 3)
560             np.array.
561 12          phase shift for each of the equivalent
562             reflections as a length (m) array
563 13          intensity correction for reflection, this times
564             :math:`F_{obs}^2` or :math:`F_{calc}^2`
565             gives Iobs or Icalc
566             (not used in single crystals?)
567 14          dict with the form factor (f or b) by atom type
568             symbol at the reflection position.
569==========  ====================================================
570
571
572*Classes and routines*
573----------------------
574
575'''
576import GSASIIpath
577GSASIIpath.SetVersionNumber("$Revision: 810 $")
578
579def LoadHistogramIDs(histList,idList):
580    '''Save the Id values for a series of histograms'''
581    VarName.IDdict['hists'] = {}
582    for h,i in zip(histList,idList):
583        VarName.IDdict['hists'][i] = h
584
585def LoadPhaseIDs(self):
586    pass
587
588class VarName(object):
589    '''Defines a GSAS-II variable either using the phase/atom/histogram
590    unique Id numbers or using a character string that specifies
591    variables by phase/atom/histogram number (which can change).
592    Note that :func:`LoadID` should be used to (re)load the current Ids
593    before creating or later using the VarName object.
594
595    A :class:`VarName` object can be created with a single parameter:
596   
597    :param str varname: a single value can be used to create a :class:`VarName`
598      object. The string must be of form "p:h:var" or "p:h:var:a", where
599
600     * p is the phase number (which may be left blank);
601     * h is the histogram number (which may be left blank);
602     * a is the atom number (which may be left blank in which case the third colon is omitted).
603
604    Alternately, a :class:`VarName` object can be created with exactly four positional parameters:
605
606    :param int phasenum: The number for the phase
607    :param int histnum: The number for the histogram
608    :param str varname: a single value can be used to create a :class:`VarName`
609    :param int atomnum: The number for the atom
610   
611    '''
612    import re
613    IDdict = {}
614    IDdict['phases'] = {}
615    IDdict['hists'] = {}
616    IDdict['atoms'] = {}
617    # This dictionary lists descriptions for GSAS-II variables.
618    # Note that keys may contain regular expressions, where '[xyz]'
619    # matches 'x' 'y' or 'z' (equivalently '[x-z]' describes this as range of values).
620    # '.*' matches any string
621    VarDesc = {
622        # Phase vars (p::<var>)
623        'A[0-5]' : 'Reciprocal metric tensor component',
624        'Vol' : 'Unit cell volume????',
625        # Atom vars (p::<var>:a)
626        'dA[xyz]' : 'change to atomic position',
627        'AUiso':'Atomic isotropic displacement parameter',
628        'AU[123][123]':'Atomic anisotropic displacement parameter',
629        'AFrac': 'Atomic occupancy parameter',
630        # Hist & Phase (HAP) vars (p:h:<var>)
631        'Bab[AU]': 'Babinet solvent scattering coef.',
632        'D[123][123]' : 'Anisotropic strain coef.',
633        'Extinction' : 'Extinction coef.',
634        'MD' : 'March-Dollase coef.',
635        'Mustrain;.*' : 'Microstrain coef.',
636        'Scale' : 'Phase scale factor',
637        'Size;.*' : 'Crystallite size value',
638        'eA' : '?',
639        #Histogram vars (:h:<var>)
640        'Absorption' : 'Absorption coef.',
641        'Displace[XY]' : 'Debye-Scherrer sample displacement',
642        'Lam' : 'Wavelength',
643        'Polariz' : 'Polarization correction',
644        'SH/L' : 'FCJ peak asymmetry correction',
645        'Scale' : 'Histogram scale factor',
646        '[UVW]' : 'Gaussian instrument broadening',
647        '[XY]' : 'Cauchy instrument broadening',
648        'Zero' : 'Debye-Scherrer zero correction',
649        'nDebye' : 'Debye model background corr. terms',
650        'nPeaks' : 'Fixed peak  background corr. terms',
651        # Global vars (::<var>)
652        }
653    def __init__(self,*args):
654        self.phase = None
655        self.histogram = None
656        self.name = ''
657        self.atom = None
658        if len(args) == 1:
659            lst = args[0].split(':')
660            raise Exception, "Need to look up IDs"
661            self.phase = lst[0]
662            self.histogram = lst[1]
663            self.name = lst[2]
664            if len(lst) > 3:
665                self.atom = lst[3]
666        elif len(args) == 4:
667            self.phase = args[0]
668            self.histogram = args[1]
669            self.name = args[2]
670            self.atom = args[3]
671        else:
672            raise Exception,"Incorrectly called GSAS-II parameter name"
673
674    def __str__(self):
675        return self.name()
676
677    def name(self):
678        '''Formats the GSAS-II variable name as a "traditional" string (p:h:<var>:a)
679
680        :returns: the variable name as a str
681        '''
682        def _fmt(val):
683            if val is None:
684                return ""
685            return str(val)
686        return _fmt(self.phase) + ":" + _fmt(self.histogram) + _fmt(self.name) + _fmt(self.atom)
687
688    def __repr__(self):
689        '''Return the detailed contents of the object
690        '''
691        s = ""
692        if self.phase:
693            s += "Phase: " + str(self.phase) + "; "
694
695        if self.histogram:
696            s += "Histogram: " + str(self.histogram) + "; "
697           
698        if self.name:
699            s += "Variable name: " + str(self.name) + "; "
700
701        if self.atom:
702            s += "Atom number: " + str(self.atom) + "; "
703
704        return s+"("+self.name()+")"
705
706    def getDescr(self):
707        '''Return a short description for a GSAS-II variable
708
709        :returns: a short description or 'no definition' if not found
710        '''
711        # iterating over uncompiled regular expressions is not terribly fast,
712        # but this routine should not need to be all that speedy
713        for key in self.VarDesc:
714            if re.match(key, self.name):
715                return self.VarDesc[key]
716        return 'no definition'
717
718    def getDescr(self):
719        '''Return a short description for a GSAS-II variable
720
721        :returns: a short description or 'no definition' if not found
722        '''
723        # iterating over uncompiled regular expressions is not terribly fast,
724        # but this routine should not need to be all that speedy
725        for key in self.VarDesc:
726            if re.match(key, self.name):
727                return self.VarDesc[key]
728        return 'no definition'
729
730    def fullDescr(self):
731        '''Return a longer description for a GSAS-II variable
732
733        :returns: a short description or 'no definition' if not found
734        '''
735        # iterating over uncompiled regular expressions is not terribly fast,
736        # but this routine should not need to be all that speedy
737        str = self.name()
738       
739        for key in self.VarDesc:
740            if re.match(key, self.name):
741                return self.VarDesc[key]
742        return 'no definition'
743
744
745    def _show(self):
746        'For testing, shows the current lookup table'
747        print 'phases', self.IDdict['phases']
748        print 'hists', self.IDdict['hists']
749        print 'atomDict', self.IDdict['atoms']
750
Note: See TracBrowser for help on using the repository browser.