source: trunk/GSASIIscriptable.py @ 3016

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

WIP

File size: 80.2 KB
Line 
1#!/usr/bin/env python
2# -*- coding: utf-8 -*-
3########### SVN repository information ###################
4# $Date: 2017-04-12 15:12:45 -0500 (Wed, 12 Apr 2017) $
5# $Author: vondreele $
6# $Revision: 2777 $
7# $URL: https://subversion.xray.aps.anl.gov/pyGSAS/trunk/GSASIIscriptable.py $
8# $Id: GSASIIIO.py 2777 2017-04-12 20:12:45Z vondreele $
9########### SVN repository information ###################
10"""
11*GSASIIscriptable: Scripting Tools*
12-----------------------------------
13
14Routines for reading, writing, modifying and creating GSAS-II project (.gpx) files.
15
16Supports a command line interface as well.
17
18Look at :class:`G2Project` to start.
19
20=====================
21Refinement parameters
22=====================
23
24There are three classes of refinement parameters:
25
26    * Histogram. Turned on and off through :func:`~G2PwdrData.set_refinements`
27      and :func:`~G2PwdrData.clear_refinements`
28    * Phase. Turned on and off through :func:`~G2Phase.set_refinements`
29      and :func:`~G2Phase.clear_refinements`
30    * Histogram-and-phase (HAP). Turned on and off through
31      :func:`~G2Phase.set_HAP_refinements` and :func:`~G2Phase.clear_HAP_refinements`
32"""
33from __future__ import division, print_function # needed?
34import argparse
35import os.path as ospath
36import datetime as dt
37import sys
38import cPickle
39import imp
40import copy
41import os
42import random as ran
43
44import numpy.ma as ma
45import scipy.interpolate as si
46import numpy as np
47import scipy as sp
48
49import GSASIIpath
50GSASIIpath.SetBinaryPath(False) # would rather have this in __name__ == '__main__' stanza
51import GSASIIobj as G2obj
52import GSASIIpwd as G2pwd
53import GSASIIstrMain as G2strMain
54import GSASIIspc as G2spc
55import GSASIIElem as G2elem
56
57# Delay imports to not slow down small scripts
58G2fil = None
59PwdrDataReaders = []
60PhaseReaders = []
61
62
63def LoadG2fil():
64    """Delay importing this module, it is slow"""
65    global G2fil
66    global PwdrDataReaders
67    global PhaseReaders
68    if G2fil is None:
69        import GSASIIfiles
70        G2fil = GSASIIfiles
71        PwdrDataReaders = G2fil.LoadImportRoutines("pwd", "Powder_Data")
72        PhaseReaders = G2fil.LoadImportRoutines("phase", "Phase")
73
74
75def LoadDictFromProjFile(ProjFile):
76    '''Read a GSAS-II project file and load items to dictionary
77    :param str ProjFile: GSAS-II project (name.gpx) full file name
78    :returns: Project,nameList, where
79
80      * Project (dict) is a representation of gpx file following the GSAS-II tree struture
81        for each item: key = tree name (e.g. 'Controls','Restraints',etc.), data is dict
82        data dict = {'data':item data whch may be list, dict or None,'subitems':subdata (if any)}
83      * nameList (list) has names of main tree entries & subentries used to reconstruct project file
84
85    Example for fap.gpx::
86
87      Project = {                 #NB:dict order is not tree order
88        u'Phases':{'data':None,'fap':{phase dict}},
89        u'PWDR FAP.XRA Bank 1':{'data':[histogram data list],'Comments':comments,'Limits':limits, etc},
90        u'Rigid bodies':{'data': {rigid body dict}},
91        u'Covariance':{'data':{covariance data dict}},
92        u'Controls':{'data':{controls data dict}},
93        u'Notebook':{'data':[notebook list]},
94        u'Restraints':{'data':{restraint data dict}},
95        u'Constraints':{'data':{constraint data dict}}]
96        }
97      nameList = [                #NB: reproduces tree order
98        [u'Notebook',],
99        [u'Controls',],
100        [u'Covariance',],
101        [u'Constraints',],
102        [u'Restraints',],
103        [u'Rigid bodies',],
104        [u'PWDR FAP.XRA Bank 1',
105             u'Comments',
106             u'Limits',
107             u'Background',
108             u'Instrument Parameters',
109             u'Sample Parameters',
110             u'Peak List',
111             u'Index Peak List',
112             u'Unit Cells List',
113             u'Reflection Lists'],
114        [u'Phases', u'fap']
115        ]
116    '''
117    # Let IOError be thrown if file does not exist
118    # if not ospath.exists(ProjFile):
119    #     print ('\n*** Error attempt to open project file that does not exist:\n   '+
120    #         str(ProjFile))
121    #     return
122    file = open(ProjFile,'rb')
123    # print('loading from file: {}'.format(ProjFile))
124    Project = {}
125    nameList = []
126    try:
127        while True:
128            try:
129                data = cPickle.load(file)
130            except EOFError:
131                break
132            datum = data[0]
133            Project[datum[0]] = {'data':datum[1]}
134            nameList.append([datum[0],])
135            for datus in data[1:]:
136                Project[datum[0]][datus[0]] = datus[1]
137                nameList[-1].append(datus[0])
138        file.close()
139        # print('project load successful')
140    except:
141        raise IOError("Error reading file "+str(ProjFile)+". This is not a GSAS-II .gpx file")
142    finally:
143        file.close()
144    return Project,nameList
145
146def SaveDictToProjFile(Project,nameList,ProjFile):
147    '''Save a GSAS-II project file from dictionary/nameList created by LoadDictFromProjFile
148
149    :param dict Project: representation of gpx file following the GSAS-II
150        tree structure as described for LoadDictFromProjFile
151    :param list nameList: names of main tree entries & subentries used to reconstruct project file
152    :param str ProjFile: full file name for output project.gpx file (including extension)
153    '''
154    file = open(ProjFile,'wb')
155    # print('save to file: {}'.format(ProjFile))
156    try:
157        for name in nameList:
158            data = []
159            item = Project[name[0]]
160            data.append([name[0],item['data']])
161            for item2 in name[1:]:
162                data.append([item2,item[item2]])
163            cPickle.dump(data,file,1)
164    finally:
165        file.close()
166    # print('project save successful')
167
168def ImportPowder(reader,filename):
169    '''Use a reader to import a powder diffraction data file
170
171    :param str reader: a scriptable reader
172    :param str filename: full name of powder data file; can be "multi-Bank" data
173
174    :returns list rdlist: list of reader objects containing powder data, one for each
175        "Bank" of data encountered in file. Items in reader object of interest are:
176
177          * rd.comments: list of str: comments found on powder file
178          * rd.dnames: list of str: data nammes suitable for use in GSASII data tree NB: duplicated in all rd entries in rdlist
179          * rd.powderdata: list of numpy arrays: pos,int,wt,zeros,zeros,zeros as needed for a PWDR entry in  GSASII data tree.
180    '''
181    rdfile,rdpath,descr = imp.find_module(reader)
182    rdclass = imp.load_module(reader,rdfile,rdpath,descr)
183    rd = rdclass.GSAS_ReaderClass()
184    if not rd.scriptable:
185        print(u'**** ERROR: '+reader+u' is not a scriptable reader')
186        return None
187    fl = open(filename,'rb')
188    rdlist = []
189    if rd.ContentsValidator(fl):
190        fl.seek(0)
191        repeat = True
192        rdbuffer = {} # create temporary storage for file reader
193        block = 0
194        while repeat: # loop if the reader asks for another pass on the file
195            block += 1
196            repeat = False
197            rd.objname = ospath.basename(filename)
198            flag = rd.Reader(filename,fl,None,buffer=rdbuffer,blocknum=block,)
199            if flag:
200                rdlist.append(copy.deepcopy(rd)) # save the result before it is written over
201                if rd.repeat:
202                    repeat = True
203        return rdlist
204    print(rd.errors)
205    return None
206
207def SetDefaultDData(dType,histoName,NShkl=0,NDij=0):
208    '''Create an initial Histogram dictionary
209
210    Author: Jackson O'Donnell (jacksonhodonnell .at. gmail.com)
211    '''
212    if dType in ['SXC','SNC']:
213        return {'Histogram':histoName,'Show':False,'Scale':[1.0,True],
214            'Babinet':{'BabA':[0.0,False],'BabU':[0.0,False]},
215            'Extinction':['Lorentzian','None', {'Tbar':0.1,'Cos2TM':0.955,
216            'Eg':[1.e-10,False],'Es':[1.e-10,False],'Ep':[1.e-10,False]}],
217            'Flack':[0.0,False]}
218    elif dType == 'SNT':
219        return {'Histogram':histoName,'Show':False,'Scale':[1.0,True],
220            'Babinet':{'BabA':[0.0,False],'BabU':[0.0,False]},
221            'Extinction':['Lorentzian','None', {
222            'Eg':[1.e-10,False],'Es':[1.e-10,False],'Ep':[1.e-10,False]}]}
223    elif 'P' in dType:
224        return {'Histogram':histoName,'Show':False,'Scale':[1.0,False],
225            'Pref.Ori.':['MD',1.0,False,[0,0,1],0,{},[],0.1],
226            'Size':['isotropic',[1.,1.,1.],[False,False,False],[0,0,1],
227                [1.,1.,1.,0.,0.,0.],6*[False,]],
228            'Mustrain':['isotropic',[1000.0,1000.0,1.0],[False,False,False],[0,0,1],
229                NShkl*[0.01,],NShkl*[False,]],
230            'HStrain':[NDij*[0.0,],NDij*[False,]],
231            'Extinction':[0.0,False],'Babinet':{'BabA':[0.0,False],'BabU':[0.0,False]}}
232
233
234def PreSetup(data):
235    '''Create part of an initial (empty) phase dictionary
236
237    from GSASIIphsGUI.py, near end of UpdatePhaseData
238
239    Author: Jackson O'Donnell (jacksonhodonnell .at. gmail.com)
240    '''
241    if 'RBModels' not in data:
242        data['RBModels'] = {}
243    if 'MCSA' not in data:
244        data['MCSA'] = {'Models':[{'Type':'MD','Coef':[1.0,False,[.8,1.2],],'axis':[0,0,1]}],'Results':[],'AtInfo':{}}
245    if 'dict' in str(type(data['MCSA']['Results'])):
246        data['MCSA']['Results'] = []
247    if 'Modulated' not in data['General']:
248        data['General']['Modulated'] = False
249    if 'modulated' in data['General']['Type']:
250        data['General']['Modulated'] = True
251        data['General']['Type'] = 'nuclear'
252
253
254def SetupGeneral(data, dirname):
255    """Helps initialize phase data.
256
257    From GSASIIphsGui.py, function of the same name. Minor changes for imports etc.
258
259    Author: Jackson O'Donnell (jacksonhodonnell .at. gmail.com)
260    """
261    mapDefault = {'MapType':'','RefList':'','Resolution':0.5,'Show bonds':True,
262                'rho':[],'rhoMax':0.,'mapSize':10.0,'cutOff':50.,'Flip':False}
263    generalData = data['General']
264    atomData = data['Atoms']
265    generalData['AtomTypes'] = []
266    generalData['Isotopes'] = {}
267
268    if 'Isotope' not in generalData:
269        generalData['Isotope'] = {}
270    if 'Data plot type' not in generalData:
271        generalData['Data plot type'] = 'Mustrain'
272    if 'POhkl' not in generalData:
273        generalData['POhkl'] = [0,0,1]
274    if 'Map' not in generalData:
275        generalData['Map'] = mapDefault.copy()
276    if 'Flip' not in generalData:
277        generalData['Flip'] = {'RefList':'','Resolution':0.5,'Norm element':'None',
278            'k-factor':0.1,'k-Max':20.,}
279    if 'testHKL' not in generalData['Flip']:
280        generalData['Flip']['testHKL'] = [[0,0,2],[2,0,0],[1,1,1],[0,2,0],[1,2,3]]
281    if 'doPawley' not in generalData:
282        generalData['doPawley'] = False     #ToDo: change to ''
283    if 'Pawley dmin' not in generalData:
284        generalData['Pawley dmin'] = 1.0
285    if 'Pawley neg wt' not in generalData:
286        generalData['Pawley neg wt'] = 0.0
287    if 'Algolrithm' in generalData.get('MCSA controls',{}) or \
288        'MCSA controls' not in generalData:
289        generalData['MCSA controls'] = {'Data source':'','Annealing':[50.,0.001,50],
290        'dmin':2.0,'Algorithm':'log','Jump coeff':[0.95,0.5],'boltzmann':1.0,
291        'fast parms':[1.0,1.0,1.0],'log slope':0.9,'Cycles':1,'Results':[],'newDmin':True}
292    if 'AtomPtrs' not in generalData:
293        generalData['AtomPtrs'] = [3,1,7,9]
294        if generalData['Type'] == 'macromolecular':
295            generalData['AtomPtrs'] = [6,4,10,12]
296        elif generalData['Type'] == 'magnetic':
297            generalData['AtomPtrs'] = [3,1,10,12]
298    if generalData['Type'] in ['modulated',]:
299        generalData['Modulated'] = True
300        generalData['Type'] = 'nuclear'
301        if 'Super' not in generalData:
302            generalData['Super'] = 1
303            generalData['SuperVec'] = [[0,0,.1],False,4]
304            generalData['SSGData'] = {}
305        if '4DmapData' not in generalData:
306            generalData['4DmapData'] = mapDefault.copy()
307            generalData['4DmapData'].update({'MapType':'Fobs'})
308    if 'Modulated' not in generalData:
309        generalData['Modulated'] = False
310    if 'HydIds' not in generalData:
311        generalData['HydIds'] = {}
312    cx,ct,cs,cia = generalData['AtomPtrs']
313    generalData['NoAtoms'] = {}
314    generalData['BondRadii'] = []
315    generalData['AngleRadii'] = []
316    generalData['vdWRadii'] = []
317    generalData['AtomMass'] = []
318    generalData['Color'] = []
319    if generalData['Type'] == 'magnetic':
320        generalData['MagDmin'] = generalData.get('MagDmin',1.0)
321        landeg = generalData.get('Lande g',[])
322    generalData['Mydir'] = dirname
323    badList = {}
324    for iat,atom in enumerate(atomData):
325        atom[ct] = atom[ct].lower().capitalize()              #force to standard form
326        if generalData['AtomTypes'].count(atom[ct]):
327            generalData['NoAtoms'][atom[ct]] += atom[cx+3]*float(atom[cs+1])
328        elif atom[ct] != 'UNK':
329            Info = G2elem.GetAtomInfo(atom[ct])
330            if not Info:
331                if atom[ct] not in badList:
332                    badList[atom[ct]] = 0
333                badList[atom[ct]] += 1
334                atom[ct] = 'UNK'
335                continue
336            atom[ct] = Info['Symbol'] # N.B. symbol might be changed by GetAtomInfo
337            generalData['AtomTypes'].append(atom[ct])
338            generalData['Z'] = Info['Z']
339            generalData['Isotopes'][atom[ct]] = Info['Isotopes']
340            generalData['BondRadii'].append(Info['Drad'])
341            generalData['AngleRadii'].append(Info['Arad'])
342            generalData['vdWRadii'].append(Info['Vdrad'])
343            if atom[ct] in generalData['Isotope']:
344                if generalData['Isotope'][atom[ct]] not in generalData['Isotopes'][atom[ct]]:
345                    isotope = generalData['Isotopes'][atom[ct]].keys()[-1]
346                    generalData['Isotope'][atom[ct]] = isotope
347                generalData['AtomMass'].append(Info['Isotopes'][generalData['Isotope'][atom[ct]]]['Mass'])
348            else:
349                generalData['Isotope'][atom[ct]] = 'Nat. Abund.'
350                if 'Nat. Abund.' not in generalData['Isotopes'][atom[ct]]:
351                    isotope = generalData['Isotopes'][atom[ct]].keys()[-1]
352                    generalData['Isotope'][atom[ct]] = isotope
353                generalData['AtomMass'].append(Info['Mass'])
354            generalData['NoAtoms'][atom[ct]] = atom[cx+3]*float(atom[cs+1])
355            generalData['Color'].append(Info['Color'])
356            if generalData['Type'] == 'magnetic':
357                if len(landeg) < len(generalData['AtomTypes']):
358                    landeg.append(2.0)
359    if generalData['Type'] == 'magnetic':
360        generalData['Lande g'] = landeg[:len(generalData['AtomTypes'])]
361
362    if badList:
363        msg = 'Warning: element symbol(s) not found:'
364        for key in badList:
365            msg += '\n\t' + key
366            if badList[key] > 1:
367                msg += ' (' + str(badList[key]) + ' times)'
368        raise Exception("Phase error:\n" + msg)
369        # wx.MessageBox(msg,caption='Element symbol error')
370    F000X = 0.
371    F000N = 0.
372    for i,elem in enumerate(generalData['AtomTypes']):
373        F000X += generalData['NoAtoms'][elem]*generalData['Z']
374        isotope = generalData['Isotope'][elem]
375        F000N += generalData['NoAtoms'][elem]*generalData['Isotopes'][elem][isotope]['SL'][0]
376    generalData['F000X'] = F000X
377    generalData['F000N'] = F000N
378    import GSASIImath as G2mth
379    generalData['Mass'] = G2mth.getMass(generalData)
380
381
382############################################
383############ CIF export helpers ############
384############################################
385## These functions are translated from
386## exports/G2export_CIF.py
387## They are used in G2Phase.export_CIF
388
389
390
391def make_empty_project(author=None, filename=None):
392    """Creates an dictionary in the style of GSASIIscriptable, for an empty
393    project.
394
395    If no author name or filename is supplied, 'no name' and
396    <current dir>/test_output.gpx are used , respectively.
397
398    Returns: project dictionary, name list
399
400    Author: Jackson O'Donnell (jacksonhodonnell .at. gmail.com)
401    """
402    if not filename:
403        filename = os.path.join(os.getcwd(), 'test_output.gpx')
404    else:
405        filename = os.path.abspath(filename)
406    gsasii_version = str(GSASIIpath.GetVersionNumber())
407    LoadG2fil()
408    import matplotlib as mpl
409    python_library_versions = G2fil.get_python_versions([mpl, np, sp])
410
411    controls_data = dict(G2obj.DefaultControls)
412    controls_data['LastSavedAs'] = unicode(filename)
413    controls_data['LastSavedUsing'] = gsasii_version
414    controls_data['PythonVersions'] = python_library_versions
415    if author:
416        controls_data['Author'] = author
417
418    output = {'Constraints': {'data': {'HAP': [], 'Hist': [], 'Phase': [],
419                                       'Global': []}},
420              'Controls': {'data': controls_data},
421              u'Covariance': {'data': {}},
422              u'Notebook': {'data': ['']},
423              u'Restraints': {'data': {}},
424              u'Rigid bodies': {'data': {'RBIds': {'Residue': [], 'Vector': []},
425                                'Residue': {'AtInfo': {}},
426                                'Vector':  {'AtInfo': {}}}}}
427
428    names = [[u'Notebook'], [u'Controls'], [u'Covariance'],
429             [u'Constraints'], [u'Restraints'], [u'Rigid bodies']]
430
431    return output, names
432
433
434class G2ImportException(Exception):
435    pass
436
437
438def import_generic(filename, readerlist):
439    """Attempt to import a filename, using a list of reader objects.
440
441    Returns the first reader object which worked."""
442    # Translated from OnImportGeneric method in GSASII.py
443    primaryReaders, secondaryReaders = [], []
444    for reader in readerlist:
445        flag = reader.ExtensionValidator(filename)
446        if flag is None:
447            secondaryReaders.append(reader)
448        elif flag:
449            primaryReaders.append(reader)
450    if not secondaryReaders and not primaryReaders:
451        raise G2ImportException("Could not read file: ", filename)
452
453    with open(filename, 'Ur') as fp:
454        rd_list = []
455
456        for rd in primaryReaders + secondaryReaders:
457            # Initialize reader
458            rd.selections = []
459            rd.dnames = []
460            rd.ReInitialize()
461            # Rewind file
462            fp.seek(0)
463            rd.errors = ""
464            if not rd.ContentsValidator(fp):
465                # Report error
466                pass
467            if len(rd.selections) > 1:
468                # Select data?
469                # GSASII.py:543
470                raise G2ImportException("Not sure what data to select")
471
472            block = 0
473            rdbuffer = {}
474            fp.seek(0)
475            repeat = True
476            while repeat:
477                repeat = False
478                block += 1
479                rd.objname = os.path.basename(filename)
480                flag = rd.Reader(filename, fp, buffer=rdbuffer, blocknum=block)
481                if flag:
482                    # Omitting image loading special cases
483                    rd.readfilename = filename
484                    rd_list.append(copy.deepcopy(rd))
485                    repeat = rd.repeat
486                else:
487                    raise G2ImportException("{}.Reader() returned:".format(rd),
488                                            flag)
489
490            if rd_list:
491                if rd.warnings:
492                    print("Read warning by", rd.formatName, "reader:",
493                          rd.warnings, file=sys.stderr)
494                return rd_list
495    raise G2ImportException("No reader could read file: " + filename)
496
497
498def load_iprms(instfile, reader):
499    """Loads instrument parameters from a file, and edits the
500    given reader.
501
502    Returns a 2-tuple of (Iparm1, Iparm2) parameters
503    """
504    LoadG2fil()
505    ext = os.path.splitext(instfile)[1]
506
507    if ext.lower() == '.instprm':
508        # New GSAS File, load appropriately
509        with open(instfile) as f:
510            lines = f.readlines()
511        bank = reader.powderentry[2]
512        numbanks = reader.numbanks
513        iparms = G2fil.ReadPowderInstprm(lines, bank, numbanks, reader)
514
515        reader.instfile = instfile
516        reader.instmsg = 'GSAS-II file' + instfile
517        return iparms
518    elif ext.lower() not in ('.prm', '.inst', '.ins'):
519        raise ValueError('Expected .prm file, found: ', instfile)
520
521    # It's an old GSAS file, load appropriately
522    Iparm = {}
523    with open(instfile, 'Ur') as fp:
524        for line in fp:
525            if '#' in line:
526                continue
527            Iparm[line[:12]] = line[12:-1]
528    ibanks = int(Iparm.get('INS   BANK  ', '1').strip())
529    if ibanks == 1:
530        reader.instbank = 1
531        reader.powderentry[2] = 1
532        reader.instfile = instfile
533        reader.instmsg = instfile + ' bank ' + str(reader.instbank)
534        return G2fil.SetPowderInstParms(Iparm, reader)
535    # TODO handle >1 banks
536    raise NotImplementedError("Check GSASIIfiles.py: ReadPowderInstprm")
537
538def load_pwd_from_reader(reader, instprm, existingnames=[]):
539    """Loads powder data from a reader object, and assembles it into a GSASII data tree.
540
541    :returns: (name, tree) - 2-tuple of the histogram name (str), and data
542
543    Author: Jackson O'Donnell (jacksonhodonnell .at. gmail.com)
544    """
545    HistName = 'PWDR ' + G2obj.StripUnicode(reader.idstring, '_')
546    HistName = unicode(G2obj.MakeUniqueLabel(HistName, existingnames))
547
548    try:
549        Iparm1, Iparm2 = instprm
550    except ValueError:
551        Iparm1, Iparm2 = load_iprms(instprm, reader)
552
553    Ymin = np.min(reader.powderdata[1])
554    Ymax = np.max(reader.powderdata[1])
555    valuesdict = {'wtFactor': 1.0,
556                  'Dummy': False,
557                  'ranId': ran.randint(0, sys.maxint),
558                  'Offset': [0.0, 0.0], 'delOffset': 0.02*Ymax,
559                  'refOffset': -0.1*Ymax, 'refDelt': 0.1*Ymax,
560                  'qPlot': False, 'dPlot': False, 'sqrtPlot': False,
561                  'Yminmax': [Ymin, Ymax]}
562    reader.Sample['ranId'] = valuesdict['ranId']
563
564    # Ending keys:
565    # [u'Reflection Lists',
566    #  u'Limits',
567    #  'data',
568    #  u'Index Peak List',
569    #  u'Comments',
570    #  u'Unit Cells List',
571    #  u'Sample Parameters',
572    #  u'Peak List',
573    #  u'Background',
574    #  u'Instrument Parameters']
575    Tmin = np.min(reader.powderdata[0])
576    Tmax = np.max(reader.powderdata[0])
577
578    default_background = [['chebyschev', False, 3, 1.0, 0.0, 0.0],
579                          {'nDebye': 0, 'debyeTerms': [], 'nPeaks': 0, 'peaksList': []}]
580
581    output_dict = {u'Reflection Lists': {},
582                   u'Limits': reader.pwdparms.get('Limits', [(Tmin, Tmax), [Tmin, Tmax]]),
583                   u'data': [valuesdict, reader.powderdata, HistName],
584                   u'Index Peak List': [[], []],
585                   u'Comments': reader.comments,
586                   u'Unit Cells List': [],
587                   u'Sample Parameters': reader.Sample,
588                   u'Peak List': {'peaks': [], 'sigDict': {}},
589                   u'Background': reader.pwdparms.get('Background', default_background),
590                   u'Instrument Parameters': [Iparm1, Iparm2],
591                   }
592
593    names = [u'Comments',
594             u'Limits',
595             u'Background',
596             u'Instrument Parameters',
597             u'Sample Parameters',
598             u'Peak List',
599             u'Index Peak List',
600             u'Unit Cells List',
601             u'Reflection Lists']
602
603    # TODO controls?? GSASII.py:1664-7
604
605    return HistName, [HistName] + names, output_dict
606
607
608def _deep_copy_into(from_, into):
609    """Helper function for reloading .gpx file. See G2Project.reload()
610
611    :author: Jackson O'Donnell (jacksonhodonnell .at. gmail.com)
612    """
613    if isinstance(from_, dict) and isinstance(into, dict):
614        combined_keys = set(from_.keys()).union(into.keys())
615        for key in combined_keys:
616            if key in from_ and key in into:
617                both_dicts = (isinstance(from_[key], dict)
618                              and isinstance(into[key], dict))
619                both_lists = (isinstance(from_[key], list)
620                              and isinstance(into[key], list))
621                if both_dicts or both_lists:
622                    _deep_copy_into(from_[key], into[key])
623                else:
624                    into[key] = from_[key]
625            elif key in from_:
626                into[key] = from_[key]
627            else:  # key in into
628                del into[key]
629    elif isinstance(from_, list) and isinstance(into, list):
630        if len(from_) == len(into):
631            for i in xrange(len(from_)):
632                both_dicts = (isinstance(from_[i], dict)
633                              and isinstance(into[i], dict))
634                both_lists = (isinstance(from_[i], list)
635                              and isinstance(into[i], list))
636                if both_dicts or both_lists:
637                    _deep_copy_into(from_[i], into[i])
638                else:
639                    into[i] = from_[i]
640        else:
641            into[:] = from_
642
643
644class G2ObjectWrapper(object):
645    """Base class for all GSAS-II object wrappers.
646
647    The underlying GSAS-II format can be accessed as `wrapper.data`. A number
648    of overrides are implemented so that the wrapper behaves like a dictionary.
649
650    Author: Jackson O'Donnell (jacksonhodonnell .at. gmail.com)
651    """
652    def __init__(self, datadict):
653        self.data = datadict
654
655    def __getitem__(self, key):
656        return self.data[key]
657
658    def __setitem__(self, key, value):
659        self.data[key] = value
660
661    def __contains__(self, key):
662        return key in self.data
663
664    def get(self, k, d=None):
665        return self.data.get(k, d)
666
667    def keys(self):
668        return self.data.keys()
669
670    def values(self):
671        return self.data.values()
672
673    def items(self):
674        return self.data.items()
675
676
677class G2Project(G2ObjectWrapper):
678    """Represents an entire GSAS-II project."""
679    def __init__(self, gpxfile=None, author=None, filename=None):
680        """Loads a GSAS-II project from a specified filename.
681
682        :param str gpxfile: Existing .gpx file to be loaded. If nonexistent,
683            creates an empty project.
684        :param str author: Author's name. Optional.
685        :param str filename: The filename the project should be saved to in
686            the future. If both filename and gpxfile are present, the project is
687            loaded from the gpxfile, then set to  save to filename in the future"""
688        if gpxfile is None:
689            filename = os.path.abspath(os.path.expanduser(filename))
690            self.filename = filename
691            self.data, self.names = make_empty_project(author=author, filename=filename)
692        elif isinstance(gpxfile, str):
693            # TODO set author, filename
694            self.filename = os.path.abspath(os.path.expanduser(gpxfile))
695            self.data, self.names = LoadDictFromProjFile(gpxfile)
696            self.index_ids()
697        else:
698            raise ValueError("Not sure what to do with gpxfile")
699
700    @classmethod
701    def from_dict_and_names(cls, gpxdict, names, filename=None):
702        """Creates a :class:`G2Project` directly from
703        a dictionary and a list of names. If in doubt, do not use this.
704
705        :returns: a :class:`G2Project`
706        """
707        out = cls()
708        if filename:
709            filename = os.path.abspath(os.path.expanduser(filename))
710            out.filename = filename
711            gpxdict['Controls']['data']['LastSavedAs'] = filename
712        else:
713            try:
714                out.filename = gpxdict['Controls']['data']['LastSavedAs']
715            except KeyError:
716                out.filename = None
717        out.data = gpxdict
718        out.names = names
719
720    def save(self, filename=None):
721        """Saves the project, either to the current filename, or to a new file.
722
723        Updates self.filename if a new filename provided"""
724        # TODO update LastSavedUsing ?
725        if filename:
726            filename = os.path.abspath(os.path.expanduser(filename))
727            self.data['Controls']['data']['LastSavedAs'] = filename
728            self.filename = filename
729        elif not self.filename:
730            raise AttributeError("No file name to save to")
731        SaveDictToProjFile(self.data, self.names, self.filename)
732
733    def add_powder_histogram(self, datafile, iparams):
734        """Loads a powder data histogram into the project.
735
736        Automatically checks for an instrument parameter file, or one can be
737        provided.
738
739        :param str datafile: The powder data file to read, a filename.
740        :param str iparams: The instrument parameters file, a filename.
741
742        :returns: A :class:`G2PwdrData` object representing
743            the histogram
744        """
745        LoadG2fil()
746        datafile = os.path.abspath(os.path.expanduser(datafile))
747        iparams = os.path.abspath(os.path.expanduser(iparams))
748        pwdrreaders = import_generic(datafile, PwdrDataReaders)
749        histname, new_names, pwdrdata = load_pwd_from_reader(
750                                          pwdrreaders[0], iparams,
751                                          [h.name for h in self.histograms()])
752        if histname in self.data:
753            print("Warning - redefining histogram", histname)
754        else:
755            if self.names[-1][0] == 'Phases':
756                self.names.insert(-1, new_names)
757            else:
758                self.names.append(new_names)
759        self.data[histname] = pwdrdata
760        return self.histogram(histname)
761
762    def add_phase(self, phasefile, phasename=None, histograms=[]):
763        """Loads a phase into the project from a .cif file
764
765        :param str phasefile: The CIF file from which to import the phase.
766        :param str phasename: The name of the new phase, or None for the default
767        :param list histograms: The names of the histograms to associate with
768            this phase
769
770        :returns: A :class:`G2Phase` object representing the
771            new phase.
772        """
773        LoadG2fil()
774        histograms = [self.histogram(h).name for h in histograms]
775        phasefile = os.path.abspath(os.path.expanduser(phasefile))
776
777        # TODO handle multiple phases in a file
778        phasereaders = import_generic(phasefile, PhaseReaders)
779        phasereader = phasereaders[0]
780
781        phasename = phasename or phasereader.Phase['General']['Name']
782        phaseNameList = [p.name for p in self.phases()]
783        phasename = G2obj.MakeUniqueLabel(phasename, phaseNameList)
784        phasereader.Phase['General']['Name'] = phasename
785
786        if 'Phases' not in self.data:
787            self.data[u'Phases'] = { 'data': None }
788        assert phasename not in self.data['Phases'], "phase names should be unique"
789        self.data['Phases'][phasename] = phasereader.Phase
790
791        if phasereader.Constraints:
792            Constraints = self.data['Constraints']
793            for i in phasereader.Constraints:
794                if isinstance(i, dict):
795                    if '_Explain' not in Constraints:
796                        Constraints['_Explain'] = {}
797                    Constraints['_Explain'].update(i)
798                else:
799                    Constraints['Phase'].append(i)
800
801        data = self.data['Phases'][phasename]
802        generalData = data['General']
803        SGData = generalData['SGData']
804        NShkl = len(G2spc.MustrainNames(SGData))
805        NDij = len(G2spc.HStrainNames(SGData))
806        Super = generalData.get('Super', 0)
807        if Super:
808            SuperVec = np.array(generalData['SuperVec'][0])
809        else:
810            SuperVec = []
811        UseList = data['Histograms']
812
813        for histname in histograms:
814            if histname.startswith('HKLF '):
815                raise NotImplementedError("Does not support HKLF yet")
816            elif histname.startswith('PWDR '):
817                hist = self.histogram(histname)
818                hist['Reflection Lists'][generalData['Name']] = {}
819                UseList[histname] = SetDefaultDData('PWDR', histname, NShkl=NShkl, NDij=NDij)
820                for key, val in [('Use', True), ('LeBail', False),
821                                 ('newLeBail', True),
822                                 ('Babinet', {'BabA': [0.0, False],
823                                              'BabU': [0.0, False]})]:
824                    if key not in UseList[histname]:
825                        UseList[histname][key] = val
826            else:
827                raise NotImplementedError("Unexpected histogram" + histname)
828
829        for obj in self.names:
830            if obj[0] == 'Phases':
831                phasenames = obj
832                break
833        else:
834            phasenames = [u'Phases']
835            self.names.append(phasenames)
836        phasenames.append(unicode(phasename))
837
838        # TODO should it be self.filename, not phasefile?
839        SetupGeneral(data, os.path.dirname(phasefile))
840        self.index_ids()
841
842        return self.phase(phasename)
843
844    def reload(self):
845        """Reload self from self.filename"""
846        data, names = LoadDictFromProjFile(self.filename)
847        self.names = names
848        # Need to deep copy the new data file data into the current tree,
849        # so that any existing G2Phase, or G2PwdrData objects will still be
850        # valid
851        _deep_copy_into(from_=data, into=self.data)
852
853    def refine(self, newfile=None, printFile=None, makeBack=False):
854        # index_ids will automatically save the project
855        self.index_ids()
856        # TODO G2strMain does not properly use printFile
857        G2strMain.Refine(self.filename, makeBack=makeBack)
858        # Reload yourself
859        self.reload()
860
861    def histogram(self, histname):
862        """Returns the histogram named histname, or None if it does not exist.
863
864        :param histname: The name of the histogram (str), or ranId or index.
865        :returns: A :class:`G2PwdrData` object, or None if
866            the histogram does not exist
867
868        .. seealso::
869            :func:`~GSASIIscriptable.G2Project.histograms`
870            :func:`~GSASIIscriptable.G2Project.phase`
871            :func:`~GSASIIscriptable.G2Project.phases`
872            """
873        if isinstance(histname, G2PwdrData):
874            if histname.proj == self:
875                return histname
876        if histname in self.data:
877            return G2PwdrData(self.data[histname], self)
878        for key, val in G2obj.HistIdLookup.items():
879            name, ranId = val
880            # histname can be either ranId (key) or index (val)
881            if ranId == histname or key == str(histname):
882                return self.histogram(name)
883
884    def histograms(self):
885        """Return a list of all histograms, as
886        :class:`G2PwdrData` objects
887
888        .. seealso::
889            :func:`~GSASIIscriptable.G2Project.histograms`
890            :func:`~GSASIIscriptable.G2Project.phase`
891            :func:`~GSASIIscriptable.G2Project.phases`
892            """
893        output = []
894        for obj in self.names:
895            if len(obj) > 1 and obj[0] != u'Phases':
896                output.append(self.histogram(obj[0]))
897        return output
898
899    def phase(self, phasename):
900        """
901        Gives an object representing the specified phase in this project.
902
903        :param str phasename: The name of the desired phase. Either the name
904            (str), the phase's ranId, or the phase's index
905        :returns: A :class:`G2Phase` object
906        :raises: KeyError
907
908        .. seealso::
909            :func:`~GSASIIscriptable.G2Project.histograms`
910            :func:`~GSASIIscriptable.G2Project.phase`
911            :func:`~GSASIIscriptable.G2Project.phases`
912            """
913        phases = self.data['Phases']
914        if phasename in phases:
915            return G2Phase(phases[phasename], phasename, self)
916        for key, val in G2obj.PhaseIdLookup.items():
917            name, ranId = val
918            # phasename can be either ranId (key) or index (val)
919            if ranId == phasename or key == str(phasename):
920                return self.phase(name)
921
922    def phases(self):
923        """
924        Returns a list of all the phases in the project.
925
926        :returns: A list of :class:`G2Phase` objects
927
928        .. seealso::
929            :func:`~GSASIIscriptable.G2Project.histogram`
930            :func:`~GSASIIscriptable.G2Project.histograms`
931            :func:`~GSASIIscriptable.G2Project.phase`
932            """
933        for obj in self.names:
934            if obj[0] == 'Phases':
935                return [self.phase(p) for p in obj[1:]]
936        return []
937
938    def do_refinements(self, refinements, histogram='all', phase='all',
939                       outputnames=None, makeBack=False):
940        """Conducts a series of refinements. Wrapper around iter_refinements
941
942        :param list refinements: A list of dictionaries defining refinements
943        :param str histogram: Name of histogram for refinements to be applied
944            to, or 'all'
945        :param str phase: Name of phase for refinements to be applied to, or
946            'all'
947        """
948        for proj in self.iter_refinements(refinements, histogram, phase,
949                                          outputnames, makeBack):
950            pass
951        return self
952
953    def iter_refinements(self, refinements, histogram='all', phase='all',
954                         outputnames=None, makeBack=False):
955        """Conducts a series of refinements, iteratively. Stops after every
956        refinement and yields this project, to allow error checking or
957        logging of intermediate results.
958
959        >>> def checked_refinements(proj):
960        ...     for p in proj.iter_refinements(refs):
961        ...         # Track intermediate results
962        ...         log(p.histogram('0').residuals)
963        ...         log(p.phase('0').get_cell())
964        ...         # Check if parameter diverged, nonsense answer, or whatever
965        ...         if is_something_wrong(p):
966        ...             raise Exception("I need a human!")
967
968        :param list refinements: A list of dictionaries defining refinements
969        :param str histogram: Name of histogram for refinements to be applied
970            to, or 'all'
971        :param str phase: Name of phase for refinements to be applied to, or
972            'all'
973        """
974        if outputnames:
975            if len(refinements) != len(outputnames):
976                raise ValueError("Should have same number of outuputs to"
977                                 "refinements")
978        else:
979            outputnames = [None for r in refinements]
980
981        for output, refinement in zip(outputnames, refinements):
982            self.set_refinement(refinement, histogram)
983            # Handle 'once' args - refinements that are disabled after this
984            # refinement
985            if 'once' in refinement:
986                temp = {'set': refinement['once']}
987                self.set_refinement(temp, histogram, phase)
988
989            if output:
990                self.save(output)
991
992            self.refine(makeBack=makeBack)
993            yield self
994
995            # Handle 'once' args - refinements that are disabled after this
996            # refinement
997            if 'once' in refinement:
998                temp = {'clear': refinement['once']}
999                self.set_refinement(temp, histogram, phase)
1000
1001    def set_refinement(self, refinement, histogram='all', phase='all'):
1002        """Apply specified refinements to a given histogram(s) or phase(s).
1003
1004        Refinement parameters are categorize in three groups:
1005
1006        1. Histogram parameters
1007        2. Phase parameters
1008        3. Histogram-and-Phase (HAP) parameters
1009
1010        :param dict refinement: The refinements to be conducted
1011        :param histogram: Either a name of a histogram (str), a list of
1012            histogram names, or 'all' (default)
1013        :param phase: Either a name of a phase (str), a list of phase names, or
1014            'all' (default)
1015
1016        .. seealso::
1017            :func:`~G2PwdrData.set_refinements`
1018            :func:`~G2PwdrData.clear_refinements`
1019            :func:`~G2Phase.set_refinements`
1020            :func:`~G2Phase.clear_refinements`
1021            :func:`~G2Phase.set_HAP_refinements`
1022            :func:`~G2Phase.clear_HAP_refinements`"""
1023
1024        if histogram == 'all':
1025            hists = self.histograms()
1026        elif isinstance(histogram, str) or isinstance(histogram, int):
1027            hists = [self.histogram(histogram)]
1028        else:
1029            hists = [self.histogram(name) for name in histogram]
1030
1031        if phase == 'all':
1032            phases = self.phases()
1033        elif isinstance(phase, str) or isinstance(phase, int):
1034            phases = [self.phase(phase)]
1035        else:
1036            phases = [self.phase(name) for name in phase]
1037
1038
1039        # TODO: HAP parameters:
1040        #   Babinet
1041        #   Extinction
1042        #   HStrain
1043        #   Mustrain
1044        #   Pref. Ori
1045        #   Size
1046
1047        pwdr_set = {}
1048        phase_set = {}
1049        hap_set = {}
1050        for key, val in refinement.get('set', {}).items():
1051            # Apply refinement options
1052            if G2PwdrData.is_valid_refinement_key(key):
1053                pwdr_set[key] = val
1054            elif G2Phase.is_valid_refinement_key(key):
1055                phase_set[key] = val
1056            elif G2Phase.is_valid_HAP_refinement_key(key):
1057                hap_set[key] = val
1058            else:
1059                raise ValueError("Unknown refinement key", key)
1060
1061        for hist in hists:
1062            hist.set_refinements(pwdr_set)
1063        for phase in phases:
1064            phase.set_refinements(phase_set)
1065        for phase in phases:
1066            phase.set_HAP_refinements(hap_set, hists)
1067
1068        pwdr_clear = {}
1069        phase_clear = {}
1070        hap_clear = {}
1071        for key, val in refinement.get('clear', {}).items():
1072            # Clear refinement options
1073            if G2PwdrData.is_valid_refinement_key(key):
1074                pwdr_clear[key] = val
1075            elif G2Phase.is_valid_refinement_key(key):
1076                phase_clear[key] = val
1077            elif G2Phase.is_valid_HAP_refinement_key(key):
1078                hap_set[key] = val
1079            else:
1080                raise ValueError("Unknown refinement key", key)
1081
1082        for hist in hists:
1083            hist.clear_refinements(pwdr_clear)
1084        for phase in phases:
1085            phase.clear_refinements(phase_clear)
1086        for phase in phases:
1087            phase.clear_HAP_refinements(hap_clear, hists)
1088
1089    def index_ids(self):
1090        import GSASIIstrIO as G2strIO
1091        self.save()
1092        return G2strIO.GetUsedHistogramsAndPhases(self.filename)
1093
1094    def add_constraint_raw(self, cons_scope, constr):
1095        """Adds a constraint of type consType to the project.
1096        cons_scope should be one of "Hist", "Phase", "HAP", or "Global".
1097
1098        WARNING it does not check the constraint is well-constructed"""
1099        constrs = self.data['Constraints']['data']
1100        if 'Global' not in constrs:
1101            constrs['Global'] = []
1102        constrs[cons_scope].append(constr)
1103
1104    def hold_many(self, vars, type):
1105        """Apply holds for all the variables in vars, for constraint of a given type.
1106
1107        type is passed directly to add_constraint_raw as consType
1108
1109        :param list vars: A list of variables to hold. Either :class:`GSASIIobj.G2VarObj` objects,
1110            string variable specifiers, or arguments for :meth:`make_var_obj`
1111        :param str type: A string constraint type specifier. See
1112        :class:`~GSASIIscriptable.G2Project.add_constraint_raw`"""
1113        for var in vars:
1114            if isinstance(var, str):
1115                var = self.make_var_obj(var)
1116            elif not isinstance(var, G2obj.G2VarObj):
1117                var = self.make_var_obj(*var)
1118            self.add_constraint_raw(type, [[1.0, var], None, None, 'h'])
1119
1120    def make_var_obj(self, phase=None, hist=None, varname=None, atomId=None,
1121                     reloadIdx=True):
1122        """Wrapper to create a G2VarObj. Takes either a string representaiton ("p:h:name:a")
1123        or individual names of phase, histogram, varname, and atomId.
1124
1125        Automatically converts string phase, hist, or atom names into the ID required
1126        by G2VarObj."""
1127
1128        if reloadIdx:
1129            self.index_ids()
1130
1131        # If string representation, short circuit
1132        if hist is None and varname is None and atomId is None:
1133            if isinstance(phase, str) and ':' in phase:
1134                return G2obj.G2VarObj(phase)
1135
1136        # Get phase index
1137        phaseObj = None
1138        if isinstance(phase, G2Phase):
1139            phaseObj = phase
1140            phase = G2obj.PhaseRanIdLookup[phase.ranId]
1141        elif phase in self.data['Phases']:
1142            phaseObj = self.phase(phase)
1143            phaseRanId = phaseObj.ranId
1144            phase = G2obj.PhaseRanIdLookup[phaseRanId]
1145
1146        # Get histogram index
1147        if isinstance(hist, G2PwdrData):
1148            hist = G2obj.HistRanIdLookup[hist.ranId]
1149        elif hist in self.data:
1150            histRanId = self.histogram(hist).ranId
1151            hist = G2obj.HistRanIdLookup[histRanId]
1152
1153        # Get atom index (if any)
1154        if isinstance(atomId, G2AtomRecord):
1155            atomId = G2obj.AtomRanIdLookup[phase][atomId.ranId]
1156        elif phaseObj:
1157            atomObj = phaseObj.atom(atomId)
1158            if atomObj:
1159                atomRanId = atomObj.ranId
1160                atomId = G2obj.AtomRanIdLookup[phase][atomRanId]
1161
1162        return G2obj.G2VarObj(phase, hist, varname, atomId)
1163
1164
1165class G2AtomRecord(G2ObjectWrapper):
1166    """Wrapper for an atom record. Has convenient accessors via @property.
1167
1168
1169    Available accessors: label, type, refinement_flags, coordinates,
1170        occupancy, ranId, id, adp_flag, uiso
1171
1172    Example:
1173
1174    >>> atom = some_phase.atom("O3")
1175    >>> # We can access the underlying data format
1176    >>> atom.data
1177    ['O3', 'O-2', '', ... ]
1178    >>> # We can also use wrapper accessors
1179    >>> atom.coordinates
1180    (0.33, 0.15, 0.5)
1181    >>> atom.refinement_flags
1182    u'FX'
1183    >>> atom.ranId
1184    4615973324315876477
1185    >>> atom.occupancy
1186    1.0
1187    """
1188    def __init__(self, data, indices, proj):
1189        self.data = data
1190        self.cx, self.ct, self.cs, self.cia = indices
1191        self.proj = proj
1192
1193    @property
1194    def label(self):
1195        return self.data[self.ct-1]
1196
1197    @property
1198    def type(self):
1199        return self.data[self.ct]
1200
1201    @property
1202    def refinement_flags(self):
1203        return self.data[self.ct+1]
1204
1205    @refinement_flags.setter
1206    def refinement_flags(self, other):
1207        # Automatically check it is a valid refinement
1208        for c in other:
1209            if c not in ' FXU':
1210                raise ValueError("Invalid atom refinement: ", other)
1211        self.data[self.ct+1] = unicode(other)
1212
1213    @property
1214    def coordinates(self):
1215        return tuple(self.data[self.cx:self.cx+3])
1216
1217    @property
1218    def occupancy(self):
1219        return self.data[self.cx+3]
1220
1221    @occupancy.setter
1222    def occupancy(self, val):
1223        self.data[self.cx+3] = float(val)
1224
1225    @property
1226    def ranId(self):
1227        return self.data[self.cia+8]
1228
1229    @property
1230    def adp_flag(self):
1231        # Either 'I' or 'A'
1232        return self.data[self.cia]
1233
1234    @property
1235    def uiso(self):
1236        if self.adp_flag == 'I':
1237            return self.data[self.cia+1]
1238        else:
1239            return self.data[self.cia+2:self.cia+8]
1240
1241
1242class G2PwdrData(G2ObjectWrapper):
1243    """Wraps a Powder Data Histogram."""
1244    def __init__(self, data, proj):
1245        self.data = data
1246        self.proj = proj
1247
1248    @staticmethod
1249    def is_valid_refinement_key(key):
1250        valid_keys = ['Limits', 'Sample Parameters', 'Background',
1251                      'Instrument Parameters']
1252        return key in valid_keys
1253
1254    @property
1255    def name(self):
1256        return self.data['data'][-1]
1257
1258    @property
1259    def ranId(self):
1260        return self.data['data'][0]['ranId']
1261
1262    @property
1263    def residuals(self):
1264        data = self.data['data'][0]
1265        return {key: data[key]
1266                for key in ['R', 'Rb', 'wR', 'wRb', 'wRmin']}
1267
1268    @property
1269    def id(self):
1270        return self.data['data'][0]['hId']
1271
1272    def fit_fixed_points(self):
1273        """Attempts to apply a background fit to the fixed points currently specified."""
1274        def SetInstParms(Inst):
1275            dataType = Inst['Type'][0]
1276            insVary = []
1277            insNames = []
1278            insVals = []
1279            for parm in Inst:
1280                insNames.append(parm)
1281                insVals.append(Inst[parm][1])
1282                if parm in ['U','V','W','X','Y','SH/L','I(L2)/I(L1)','alpha',
1283                    'beta-0','beta-1','beta-q','sig-0','sig-1','sig-2','sig-q',] and Inst[parm][2]:
1284                        Inst[parm][2] = False
1285            instDict = dict(zip(insNames, insVals))
1286            instDict['X'] = max(instDict['X'], 0.01)
1287            instDict['Y'] = max(instDict['Y'], 0.01)
1288            if 'SH/L' in instDict:
1289                instDict['SH/L'] = max(instDict['SH/L'], 0.002)
1290            return dataType, instDict, insVary
1291
1292        bgrnd = self.data['Background']
1293
1294        # Need our fixed points in order
1295        bgrnd[1]['FixedPoints'].sort(key=lambda pair: pair[0])
1296        X = [x for x, y in bgrnd[1]['FixedPoints']]
1297        Y = [y for x, y in bgrnd[1]['FixedPoints']]
1298
1299        limits = self.data['Limits'][1]
1300        if X[0] > limits[0]:
1301            X = [limits[0]] + X
1302            Y = [Y[0]] + Y
1303        if X[-1] < limits[1]:
1304            X += [limits[1]]
1305            Y += [Y[-1]]
1306
1307        # Some simple lookups
1308        controls = self.proj['Controls']['data']
1309        inst, inst2 = self.data['Instrument Parameters']
1310        pwddata = self.data['data'][1]
1311
1312        # Construct the data for background fitting
1313        xBeg = np.searchsorted(pwddata[0], limits[0])
1314        xFin = np.searchsorted(pwddata[0], limits[1])
1315        xdata = pwddata[0][xBeg:xFin]
1316        ydata = si.interp1d(X,Y)(ma.getdata(xdata))
1317
1318        W = [1]*len(xdata)
1319        Z = [0]*len(xdata)
1320
1321        dataType, insDict, insVary = SetInstParms(inst)
1322        bakType, bakDict, bakVary = G2pwd.SetBackgroundParms(bgrnd)
1323
1324        # Do the fit
1325        data = np.array([xdata, ydata, W, Z, Z, Z])
1326        G2pwd.DoPeakFit('LSQ', [], bgrnd, limits, inst, inst2, data,
1327                        prevVaryList=bakVary, controls=controls)
1328
1329        # Post-fit
1330        parmDict = {}
1331        bakType, bakDict, bakVary = G2pwd.SetBackgroundParms(bgrnd)
1332        parmDict.update(bakDict)
1333        parmDict.update(insDict)
1334        pwddata[3][xBeg:xFin] *= 0
1335        pwddata[5][xBeg:xFin] *= 0
1336        pwddata[4][xBeg:xFin] = G2pwd.getBackground('', parmDict, bakType, dataType, xdata)[0]
1337
1338        # TODO adjust pwddata? GSASIIpwdGUI.py:1041
1339        # TODO update background
1340        self.proj.save()
1341
1342    def y_calc(self):
1343        return self.data['data'][1][3]
1344
1345    def plot(self, Yobs=True, Ycalc=True, Background=True, Residual=True):
1346        try:
1347            import matplotlib.pyplot as plt
1348            data = self.data['data'][1]
1349            if Yobs:
1350                plt.plot(data[0], data[1], label='Yobs')
1351            if Ycalc:
1352                plt.plot(data[0], data[3], label='Ycalc')
1353            if Background:
1354                plt.plot(data[0], data[4], label='Background')
1355            if Residual:
1356                plt.plot(data[0], data[5], label="Residual")
1357        except ImportError:
1358            pass
1359
1360    def set_refinements(self, refs):
1361        """Sets the refinement parameter 'key' to the specification 'value'
1362
1363        :param dict refs: A dictionary of the parameters to be set.
1364
1365        :returns: None"""
1366        do_fit_fixed_points = False
1367        for key, value in refs.items():
1368            if key == 'Limits':
1369                old_limits = self.data['Limits'][1]
1370                new_limits = value
1371                if isinstance(new_limits, dict):
1372                    if 'low' in new_limits:
1373                        old_limits[0] = new_limits['low']
1374                    if 'high' in new_limits:
1375                        old_limits[1] = new_limits['high']
1376                else:
1377                    old_limits[0], old_limits[1] = new_limits
1378            elif key == 'Sample Parameters':
1379                sample = self.data['Sample Parameters']
1380                for sparam in value:
1381                    sample[sparam][1] = True
1382            elif key == 'Background':
1383                bkg, peaks = self.data['Background']
1384
1385                # If True or False, just set the refine parameter
1386                if value in (True, False):
1387                    bkg[1] = value
1388                    return
1389
1390                if 'type' in value:
1391                    bkg[0] = value['type']
1392                if 'refine' in value:
1393                    bkg[1] = value['refine']
1394                if 'no. coeffs' in value:
1395                    cur_coeffs = bkg[2]
1396                    n_coeffs = value['no. coeffs']
1397                    if n_coeffs > cur_coeffs:
1398                        for x in range(n_coeffs - cur_coeffs):
1399                            bkg.append(0.0)
1400                    else:
1401                        for _ in range(cur_coeffs - n_coeffs):
1402                            bkg.pop()
1403                    bkg[2] = n_coeffs
1404                if 'coeffs' in value:
1405                    bkg[3:] = value['coeffs']
1406                if 'FixedPoints' in value:
1407                    peaks['FixedPoints'] = [(float(a), float(b))
1408                                            for a, b in value['FixedPoints']]
1409                if value.get('fit fixed points', False):
1410                    do_fit_fixed_points = True
1411
1412            elif key == 'Instrument Parameters':
1413                instrument, secondary = self.data['Instrument Parameters']
1414                for iparam in value:
1415                    try:
1416                        instrument[iparam][2] = True
1417                    except IndexError:
1418                        raise ValueError("Invalid key:", iparam)
1419            else:
1420                raise ValueError("Unknown key:", key)
1421        # Fit fixed points after the fact - ensure they are after fixed points
1422        # are added
1423        if do_fit_fixed_points:
1424            # Background won't be fit if refinement flag not set
1425            orig = self.data['Background'][0][1]
1426            self.data['Background'][0][1] = True
1427            self.fit_fixed_points()
1428            # Restore the previous value
1429            self.data['Background'][0][1] = orig
1430
1431    def clear_refinements(self, refs):
1432        """Clears the refinement parameter 'key' and its associated value.
1433
1434        :param dict refs: A dictionary of parameters to clear."""
1435        for key, value in refs.items():
1436            if key == 'Limits':
1437                old_limits, cur_limits = self.data['Limits']
1438                cur_limits[0], cur_limits[1] = old_limits
1439            elif key == 'Sample Parameters':
1440                sample = self.data['Sample Parameters']
1441                for sparam in value:
1442                    sample[sparam][1] = False
1443            elif key == 'Background':
1444                bkg, peaks = self.data['Background']
1445
1446                # If True or False, just set the refine parameter
1447                if value in (True, False):
1448                    bkg[1] = False
1449                    return
1450
1451                bkg[1] = False
1452                if 'FixedPoints' in value:
1453                    if 'FixedPoints' in peaks:
1454                        del peaks['FixedPoints']
1455            elif key == 'Instrument Parameters':
1456                instrument, secondary = self.data['Instrument Parameters']
1457                for iparam in value:
1458                    instrument[iparam][2] = False
1459            else:
1460                raise ValueError("Unknown key:", key)
1461
1462
1463class G2Phase(G2ObjectWrapper):
1464    """A wrapper object around a given phase.
1465
1466    Author: Jackson O'Donnell (jacksonhodonnell .at. gmail.com)
1467    """
1468    def __init__(self, data, name, proj):
1469        self.data = data
1470        self.name = name
1471        self.proj = proj
1472
1473    @staticmethod
1474    def is_valid_refinement_key(key):
1475        valid_keys = ["Cell", "Atoms", "LeBail"]
1476        return key in valid_keys
1477
1478    @staticmethod
1479    def is_valid_HAP_refinement_key(key):
1480        valid_keys = ["Babinet", "Extinction", "HStrain", "Mustrain",
1481                      "Pref.Ori.", "Show", "Size", "Use", "Scale"]
1482        return key in valid_keys
1483
1484    def atom(self, atomlabel):
1485        """Returns the atom specified by atomlabel, or None if it does not
1486        exist.
1487
1488        :param str atomlabel: The name of the atom (e.g. "O2")
1489        :returns: A :class:`G2AtomRecord` object
1490            representing the atom.
1491        """
1492        # Consult GSASIIobj.py for the meaning of this
1493        cx, ct, cs, cia = self.data['General']['AtomPtrs']
1494        ptrs = [cx, ct, cs, cia]
1495        atoms = self.data['Atoms']
1496        for atom in atoms:
1497            if atom[ct-1] == atomlabel:
1498                return G2AtomRecord(atom, ptrs, self.proj)
1499
1500    def atoms(self):
1501        """Returns a list of atoms present in the phase.
1502
1503        :returns: A list of :class:`G2AtomRecord` objects.
1504
1505        .. seealso::
1506            :func:`~GSASIIscriptable.G2Phase.atom`
1507            :class:`G2AtomRecord`
1508        """
1509        ptrs = self.data['General']['AtomPtrs']
1510        output = []
1511        atoms = self.data['Atoms']
1512        for atom in atoms:
1513            output.append(G2AtomRecord(atom, ptrs, self.proj))
1514        return output
1515
1516    def histograms(self):
1517        output = []
1518        for hname in self.data.get('Histograms', {}).keys():
1519            output.append(self.proj.histogram(hname))
1520        return output
1521
1522    @property
1523    def ranId(self):
1524        return self.data['ranId']
1525
1526    @property
1527    def id(self):
1528        return self.data['pId']
1529
1530    def get_cell(self):
1531        """Returns a dictionary of the cell parameters, with keys:
1532            'length_a', 'length_b', 'length_c', 'angle_alpha', 'angle_beta', 'angle_gamma', 'volume'
1533
1534        :returns: a dict
1535
1536        .. seealso::
1537           :func:`~G2Phase.get_cell_and_esd`
1538        """
1539        cell = self.data['General']['Cell']
1540        return {'length_a': cell[1], 'length_b': cell[2], 'length_c': cell[3],
1541                'angle_alpha': cell[4], 'angle_beta': cell[5], 'angle_gamma': cell[6],
1542                'volume': cell[7]}
1543
1544    def get_cell_and_esd(self):
1545        """
1546        Returns a pair of dictionaries, the first representing the unit cell, the second
1547        representing the estimated standard deviations of the unit cell.
1548
1549        :returns: a tuple of two dictionaries
1550        .. seealso::
1551           :func:`~G2Phase.get_cell`
1552        """
1553        # translated from GSASIIstrIO.ExportBaseclass.GetCell
1554        import GSASIIstrIO as G2stIO
1555        import GSASIIlattice as G2lat
1556        import GSASIImapvars as G2mv
1557        try:
1558            pfx = str(self.id) + '::'
1559            sgdata = self['General']['SGData']
1560            covDict = self.proj['Covariance']['data']
1561
1562            parmDict = dict(zip(covDict.get('varyList',[]),
1563                                covDict.get('variables',[])))
1564            sigDict = dict(zip(covDict.get('varyList',[]),
1565                               covDict.get('sig',[])))
1566
1567            if covDict.get('covMatrix') is not None:
1568                sigDict.update(G2mv.ComputeDepESD(covDict['covMatrix'],
1569                                                  covDict['varyList'],
1570                                                  parmDict))
1571
1572            A, sigA = G2stIO.cellFill(pfx, sgdata, parmDict, sigDict)
1573            cellSig = G2stIO.getCellEsd(pfx, sgdata, A, self.proj['Covariance']['data'])
1574            cellList = G2lat.A2cell(A) + (G2lat.calc_V(A),)
1575            cellDict, cellSigDict = {}, {}
1576            for i, key in enumerate(['length_a', 'length_b', 'length_c',
1577                                     'angle_alpha', 'angle_beta', 'angle_gamma',
1578                                     'volume']):
1579                cellDict[key] = cellList[i]
1580                cellSigDict[key] = cellSig[i]
1581            return cellDict, cellSigDict
1582        except KeyError:
1583            cell = self.get_cell()
1584            return cell, {key: 0.0 for key in cell}
1585
1586    def export_CIF(self, outputname, quickmode=True):
1587        """Write this phase to a .cif file named outputname
1588
1589        :param str outputname: The name of the .cif file to write to
1590        :param bool quickmode: Currently ignored. Carryover from exports.G2export_CIF"""
1591        # This code is all taken from exports/G2export_CIF.py
1592        # Functions copied have the same names
1593        import GSASIImath as G2mth
1594        import GSASIImapvars as G2mv
1595        from exports import G2export_CIF as cif
1596
1597        CIFdate = dt.datetime.strftime(dt.datetime.now(),"%Y-%m-%dT%H:%M")
1598        CIFname = os.path.splitext(self.proj.filename)[0]
1599        CIFname = os.path.split(CIFname)[1]
1600        CIFname = ''.join([c if ord(c) < 128 else ''
1601                           for c in CIFname.replace(' ', '_')])
1602        try:
1603            author = self.proj['Controls']['data'].get('Author','').strip()
1604        except KeyError:
1605            pass
1606        oneblock = True
1607
1608        covDict = self.proj['Covariance']['data']
1609        parmDict = dict(zip(covDict.get('varyList',[]),
1610                            covDict.get('variables',[])))
1611        sigDict = dict(zip(covDict.get('varyList',[]),
1612                           covDict.get('sig',[])))
1613
1614        if covDict.get('covMatrix') is not None:
1615            sigDict.update(G2mv.ComputeDepESD(covDict['covMatrix'],
1616                                              covDict['varyList'],
1617                                              parmDict))
1618
1619        with open(outputname, 'w') as fp:
1620            fp.write(' \n' + 70*'#' + '\n')
1621            cif.WriteCIFitem(fp, 'data_' + CIFname)
1622            # from exports.G2export_CIF.WritePhaseInfo
1623            cif.WriteCIFitem(fp, '\n# phase info for '+str(self.name) + ' follows')
1624            cif.WriteCIFitem(fp, '_pd_phase_name', self.name)
1625            # TODO get esds
1626            cellDict = self.get_cell()
1627            defsigL = 3*[-0.00001] + 3*[-0.001] + [-0.01] # significance to use when no sigma
1628            names = ['length_a','length_b','length_c',
1629                     'angle_alpha','angle_beta ','angle_gamma',
1630                     'volume']
1631            for key, val in cellDict.items():
1632                cif.WriteCIFitem(fp, '_cell_' + key, G2mth.ValEsd(val))
1633
1634            cif.WriteCIFitem(fp, '_symmetry_cell_setting',
1635                         self.data['General']['SGData']['SGSys'])
1636
1637            spacegroup = self.data['General']['SGData']['SpGrp'].strip()
1638            # regularize capitalization and remove trailing H/R
1639            spacegroup = spacegroup[0].upper() + spacegroup[1:].lower().rstrip('rh ')
1640            cif.WriteCIFitem(fp, '_symmetry_space_group_name_H-M', spacegroup)
1641
1642            # generate symmetry operations including centering and center of symmetry
1643            SymOpList, offsetList, symOpList, G2oprList, G2opcodes = G2spc.AllOps(
1644                self.data['General']['SGData'])
1645            cif.WriteCIFitem(fp, 'loop_\n    _space_group_symop_id\n    _space_group_symop_operation_xyz')
1646            for i, op in enumerate(SymOpList,start=1):
1647                cif.WriteCIFitem(fp, '   {:3d}  {:}'.format(i,op.lower()))
1648
1649            # TODO skipped histograms, exports/G2export_CIF.py:880
1650
1651            # report atom params
1652            if self.data['General']['Type'] in ['nuclear','macromolecular']:        #this needs macromolecular variant, etc!
1653                cif.WriteAtomsNuclear(fp, self.data, self.name, parmDict, sigDict, [])
1654                # self._WriteAtomsNuclear(fp, parmDict, sigDict)
1655            else:
1656                raise Exception,"no export for "+str(self.data['General']['Type'])+" coordinates implemented"
1657            # report cell contents
1658            cif.WriteComposition(fp, self.data, self.name, parmDict)
1659            if not quickmode and self.data['General']['Type'] == 'nuclear':      # report distances and angles
1660                # WriteDistances(fp,self.name,SymOpList,offsetList,symOpList,G2oprList)
1661                raise NotImplementedError("only quickmode currently supported")
1662            if 'Map' in self.data['General'] and 'minmax' in self.data['General']['Map']:
1663                cif.WriteCIFitem(fp,'\n# Difference density results')
1664                MinMax = self.data['General']['Map']['minmax']
1665                cif.WriteCIFitem(fp,'_refine_diff_density_max',G2mth.ValEsd(MinMax[0],-0.009))
1666                cif.WriteCIFitem(fp,'_refine_diff_density_min',G2mth.ValEsd(MinMax[1],-0.009))
1667
1668
1669    def set_refinements(self, refs):
1670        """Sets the refinement parameter 'key' to the specification 'value'
1671
1672        :param dict refs: A dictionary of the parameters to be set.
1673
1674        :returns: None"""
1675        for key, value in refs.items():
1676            if key == "Cell":
1677                self.data['General']['Cell'][0] = True
1678            elif key == "Atoms":
1679                cx, ct, cs, cia = self.data['General']['AtomPtrs']
1680
1681                for atomlabel, atomrefinement in value.items():
1682                    if atomlabel == 'all':
1683                        for atom in self.atoms():
1684                            atom.refinement_flags = atomrefinement
1685                    else:
1686                        atom = self.atom(atomlabel)
1687                        atom.refinement_flags = atomrefinement
1688            elif key == "LeBail":
1689                hists = self.data['Histograms']
1690                for hname, hoptions in hists.items():
1691                    if 'LeBail' not in hoptions:
1692                        hoptions['newLeBail'] = True
1693                    hoptions['LeBail'] = bool(value)
1694            else:
1695                raise ValueError("Unknown key:", key)
1696
1697    def clear_refinements(self, refs):
1698        """Clears a given set of parameters.
1699
1700        :param dict refs: The parameters to clear"""
1701        for key, value in refs.items():
1702            if key == "Cell":
1703                self.data['General']['Cell'][0] = False
1704            elif key == "Atoms":
1705                cx, ct, cs, cia = self.data['General']['AtomPtrs']
1706
1707                for atomlabel in value:
1708                    atom = self.atom(atomlabel)
1709                    # Set refinement to none
1710                    atom.refinement_flags = ' '
1711            elif key == "LeBail":
1712                hists = self.data['Histograms']
1713                for hname, hoptions in hists.items():
1714                    if 'LeBail' not in hoptions:
1715                        hoptions['newLeBail'] = True
1716                    hoptions['LeBail'] = False
1717            else:
1718                raise ValueError("Unknown key:", key)
1719
1720    def set_HAP_refinements(self, refs, histograms='all'):
1721        """Sets the given HAP refinement parameters between this phase and
1722        the given histograms
1723
1724        :param dict refs: A dictionary of the parameters to be set.
1725        :param histograms: Either 'all' (default) or a list of the histograms
1726            whose HAP parameters will be set with this phase. Histogram and phase
1727            must already be associated
1728
1729        :returns: None
1730        """
1731        if histograms == 'all':
1732            histograms = self.data['Histograms'].values()
1733        else:
1734            histograms = [self.data['Histograms'][h.name] for h in histograms]
1735
1736        for key, val in refs.items():
1737            for h in histograms:
1738                if key == 'Babinet':
1739                    try:
1740                        sets = list(val)
1741                    except ValueError:
1742                        sets = ['BabA', 'BabU']
1743                    for param in sets:
1744                        if param not in ['BabA', 'BabU']:
1745                            raise ValueError("Not sure what to do with" + param)
1746                        for hist in histograms:
1747                            hist['Babinet'][param][1] = True
1748                elif key == 'Extinction':
1749                    for h in histograms:
1750                        h['Extinction'][1] = True
1751                elif key == 'HStrain':
1752                    for h in histograms:
1753                        hist['HStrain'][1] = [True for p in hist['Hstrain'][0]]
1754                elif key == 'Mustrain':
1755                    for h in histograms:
1756                        mustrain = h['Mustrain']
1757                        newType = None
1758                        if isinstance(val, str):
1759                            if val in ['isotropic', 'uniaxial', 'generalized']:
1760                                newType = val
1761                            else:
1762                                raise ValueError("Not a Mustrain type: " + val)
1763                        elif isinstance(val, dict):
1764                            newType = val.get('type', None)
1765                            direction = val.get('direction', None)
1766
1767                        if newType:
1768                            mustrain[0] = newType
1769                            if newType == 'isotropic':
1770                                mustrain[2][0] = True
1771                                mustrain[5] = [False for p in mustrain[4]]
1772                            elif newType == 'uniaxial':
1773                                pass
1774                            else:  # newtype == 'generalized'
1775                                mustrain[2] = [False for p in mustrain[1]]
1776                        if direction:
1777                            direction = [int(n) for n in direction]
1778                            if len(direction) != 3:
1779                                raise ValueError("Expected hkl, found", direction)
1780                            mustrain[3] = direction
1781                elif key == 'Pref.Ori.':
1782                    for h in histograms:
1783                        h['Pref.Ori.'][2] = True
1784                elif key == 'Show':
1785                    for h in histograms:
1786                        h['Show'] = True
1787                elif key == 'Size':
1788                    raise NotImplementedError()
1789                elif key == 'Use':
1790                    for h in histograms:
1791                        h['Use'] = True
1792                elif key == 'Scale':
1793                    for h in histograms:
1794                        h['Scale'][1] = False
1795
1796    def clear_HAP_refinements(self, refs, histograms='all'):
1797        """Clears the given HAP refinement parameters between this phase and
1798        the given histograms
1799
1800        :param dict refs: A dictionary of the parameters to be cleared.
1801        :param histograms: Either 'all' (default) or a list of the histograms
1802            whose HAP parameters will be cleared with this phase. Histogram and
1803            phase must already be associated
1804
1805        :returns: None
1806        """
1807        if histograms == 'all':
1808            histograms = self.data['Histograms'].values()
1809        else:
1810            histograms = [self.data['Histograms'][h.name] for h in histograms]
1811
1812        for key, val in refs.items():
1813            for h in histograms:
1814                if key == 'Babinet':
1815                    try:
1816                        sets = list(val)
1817                    except ValueError:
1818                        sets = ['BabA', 'BabU']
1819                    for param in sets:
1820                        if param not in ['BabA', 'BabU']:
1821                            raise ValueError("Not sure what to do with" + param)
1822                        for hist in histograms:
1823                            hist['Babinet'][param][1] = False
1824                elif key == 'Extinction':
1825                    for h in histograms:
1826                        h['Extinction'][1] = False
1827                elif key == 'HStrain':
1828                    for h in histograms:
1829                        hist['HStrain'][1] = [False for p in hist['Hstrain'][0]]
1830                elif key == 'Mustrain':
1831                    for h in histograms:
1832                        mustrain = h['Mustrain']
1833                        mustrain[2] = [False for p in mustrain[1]]
1834                        mustrain[5] = [False for p in mustrain[4]]
1835                elif key == 'Pref.Ori.':
1836                    for h in histograms:
1837                        h['Pref.Ori.'][2] = False
1838                elif key == 'Show':
1839                    for h in histograms:
1840                        h['Show'] = False
1841                elif key == 'Size':
1842                    raise NotImplementedError()
1843                elif key == 'Use':
1844                    for h in histograms:
1845                        h['Use'] = False
1846                elif key == 'Scale':
1847                    for h in histograms:
1848                        h['Scale'][1] = False
1849
1850
1851##########################
1852# Command Line Interface #
1853##########################
1854
1855
1856def create(args):
1857    """The create subcommand."""
1858    proj = G2Project(filename=args.filename)
1859
1860    hist_objs = []
1861    for h in args.histograms:
1862        hist_objs.append(proj.add_powder_histogram(h, args.iparams))
1863
1864    for p in args.phases:
1865        proj.add_phase(p, histograms=hist_objs)
1866
1867    proj.save()
1868
1869
1870def dump(args):
1871    """The dump subcommand"""
1872    if not args.histograms and not args.phases:
1873        args.raw = True
1874    if args.raw:
1875        import IPython.lib.pretty as pretty
1876
1877    for fname in args.files:
1878        if args.raw:
1879            proj, nameList = LoadDictFromProjFile(fname)
1880            print("file:", fname)
1881            print("names:", nameList)
1882            for key, val in proj.items():
1883                print(key, ":")
1884                pretty.pprint(val)
1885        else:
1886            proj = G2Project(fname)
1887            if args.histograms:
1888                hists = proj.histograms()
1889                for h in hists:
1890                    print(fname, "hist", h.id, h.name)
1891            if args.phases:
1892                phase_list = proj.phases()
1893                for p in phase_list:
1894                    print(fname, "phase", p.id, p.name)
1895
1896
1897def IPyBrowse(args):
1898    """Load a .gpx file and then open a IPython shell to browse it
1899    """
1900    for fname in args.files:
1901        proj, nameList = LoadDictFromProjFile(fname)
1902        msg = "\nfname {} loaded into proj (dict) with names in nameList".format(fname)
1903        GSASIIpath.IPyBreak_base(msg)
1904        break
1905
1906
1907def refine(args):
1908    """The refine subcommand"""
1909    proj = G2Project(args.gpxfile)
1910    if args.refinements is None:
1911        proj.refine()
1912    else:
1913        import json
1914        with open(args.refinements) as refs:
1915            refs = json.load(refs)
1916        proj.do_refinements(refs['refinements'])
1917
1918
1919def seqrefine(args):
1920    """The seqrefine subcommand"""
1921    raise NotImplementedError("seqrefine is not yet implemented")
1922
1923
1924def export(args):
1925    """The export subcommand"""
1926    # Export phase as CIF to args.exportfile
1927    proj = G2Project(args.gpxfile)
1928    phase = proj.phase(args.phase)
1929    phase.export_CIF(args.exportfile)
1930
1931
1932def _args_kwargs(*args, **kwargs):
1933    return args, kwargs
1934
1935# A dictionary of the name of each subcommand, and a tuple
1936# of its associated function and a list of its arguments
1937# The arguments are passed directly to the add_argument() method
1938# of an argparse.ArgumentParser
1939subcommands = {"create":
1940               (create, [_args_kwargs('filename',
1941                                      help='the project file to create. should end in .gpx'),
1942
1943                         _args_kwargs('-i', '--iparams',
1944                                      help='instrument parameter file'),
1945
1946                         _args_kwargs('-d', '--histograms',
1947                                      nargs='+',
1948                                      help='list of datafiles to add as histograms'),
1949
1950                         _args_kwargs('-p', '--phases',
1951                                      nargs='+',
1952                                      help='list of phases to add. phases are '
1953                                           'automatically associated with all '
1954                                           'histograms given.')]),
1955
1956               "dump": (dump, [_args_kwargs('-d', '--histograms',
1957                                     action='store_true',
1958                                     help='list histograms in files, overrides --raw'),
1959
1960                               _args_kwargs('-p', '--phases',
1961                                            action='store_true',
1962                                            help='list phases in files, overrides --raw'),
1963
1964                               _args_kwargs('-r', '--raw',
1965                                      action='store_true', help='dump raw file contents, default'),
1966
1967                               _args_kwargs('files', nargs='+')]),
1968
1969               "refine":
1970               (refine, [_args_kwargs('gpxfile', help='the project file to refine'),
1971                         _args_kwargs('refinements',
1972                                      help='json file of refinements to apply. if not present'
1973                                           ' refines file as-is',
1974                                      default=None,
1975                                      nargs='?')]),
1976
1977               "seqrefine": (seqrefine, []),
1978               "export": (export, [_args_kwargs('gpxfile',
1979                                                help='the project file from which to export'),
1980                                   _args_kwargs('phase', help='identifier of phase to export'),
1981                                   _args_kwargs('exportfile', help='the .cif file to export to')]),
1982               "browse": (IPyBrowse, [_args_kwargs('files', nargs='+',
1983                                                   help='list of files to browse')])}
1984
1985
1986def main():
1987    '''The command line interface for GSASIIscriptable.
1988
1989    Executes one of the following subcommands:
1990
1991        * :func:`create`
1992        * :func:`dump`
1993        * :func:`refine`
1994        * :func:`seqrefine`
1995        * :func:`export`
1996        * :func:`browse`
1997
1998    .. seealso::
1999        :func:`create`
2000        :func:`dump`
2001        :func:`refine`
2002        :func:`seqrefine`
2003        :func:`export`
2004        :func:`browse`
2005    '''
2006    parser = argparse.ArgumentParser()
2007    subs = parser.add_subparsers()
2008
2009    # Create all of the specified subparsers
2010    for name, (func, args) in subcommands.items():
2011        new_parser = subs.add_parser(name)
2012        for listargs, kwds in args:
2013            new_parser.add_argument(*listargs, **kwds)
2014        new_parser.set_defaults(func=func)
2015
2016    # Parse and trigger subcommand
2017    result = parser.parse_args()
2018    result.func(result)
2019
2020    # argv = sys.argv
2021    # if len(argv) > 1 and argv[1] in subcommands:
2022    #     subcommands[argv[1]](*argv[2:])
2023    # elif len(argv) == 1 or argv[1] in ('help', '--help', '-h'):
2024    #     # TODO print usage
2025    #     subcommand_names = ' | '.join(sorted(subcommands.keys()))
2026    #     print("USAGE: {} [ {} ] ...".format(argv[0], subcommand_names))
2027    # else:
2028    #     print("Unknown subcommand: {}".format(argv[1]))
2029    #     print("Available subcommands:")
2030    #     for name in sorted(subcommands.keys()):
2031    #         print("\t{}".format(name))
2032    #     sys.exit(-1)
2033    # sys.exit(0)
2034
2035if __name__ == '__main__':
2036    main()
2037
2038
2039    # from refinements.py
2040#      USAGE = """USAGE: {} datafile instparams phasefile projectname refinements
2041
2042# datafile:    Input powder data
2043# intparams:   Corresponding instrument parameter file
2044# phasefile:   Phase to refine against data
2045# projectname: Project file to be created, should end in .gpx
2046# refinements: JSON file of refinements to be executed
2047# """
2048#     try:
2049#         datafile, instprm, phasefile, projectname, refinements = sys.argv[1:]
2050#     except ValueError:
2051#         print(USAGE.format(sys.argv[0]))
2052#         sys.exit(-1)
2053
2054#     try:
2055#         with open(refinements) as f:
2056#             refinements = json.load(f)
2057#     except IOError:
2058#         print("No such refinements file: {}".format(refinements))
2059
2060#     print("Creating project file \"{}\"...".format(projectname))
2061#     proj = G2Project(filename=projectname)
2062#     # Add the histogram
2063#     hist = proj.add_powder_histogram(datafile, instprm)
2064#     # Add the phase, and associate it with the histogram
2065#     proj.add_phase(phasefile, histograms=[hist.name])
2066
2067#     proj.do_refinements(refinements['refinements'])
2068#     proj.save()
2069
2070
2071    # from gpx_dumper
2072    # import IPython.lib.pretty as pretty
2073    # proj, nameList = LoadDictFromProjFile(sys.argv[1])
2074    # print("names:", nameList)
2075    # for key, val in proj.items():
2076    #     print(key, ":", sep='')
2077    #     pretty.pprint(val)
Note: See TracBrowser for help on using the repository browser.