Changeset 2963 for branch


Ignore:
Timestamp:
Aug 6, 2017 3:08:30 PM (6 years ago)
Author:
toby
Message:

Incorporate Jackson O'Donnell's gpx_manipulations code

Location:
branch/2frame
Files:
1 added
4 edited

Legend:

Unmodified
Added
Removed

  • branch/2frame/GSASIIIO.py

    r2915 r2963  
    7979    '''
    8080    return re.sub('\s+', ' ', val).strip()
    81 
    82 def makeInstDict(names,data,codes):
    83     inst = dict(zip(names,zip(data,data,codes)))
    84     for item in inst:
    85         inst[item] = list(inst[item])
    86     return inst
    8781
    8882def FileDlgFixExt(dlg,file):

  • branch/2frame/GSASIIdataGUI.py

    r2957 r2963  
    4040import GSASIImath as G2mth
    4141import GSASIIIO as G2IO
     42import GSASIIfiles as G2fil
    4243import GSASIIstrIO as G2stIO
    4344import GSASIIlattice as G2lat
    … …  
    392393        self.SetTopWindow(self.main)
    393394        # save the current package versions
    394         self.main.PackageVersions = []
    395         self.main.PackageVersions.append(['Python',sys.version.split()[0]])
    396         for p in (wx,mpl,np,sp,ogl):
    397             self.main.PackageVersions.append([p.__name__,p.__version__])
    398         # TODO: not sure how this works with pillow
    399         try:
    400             from PIL import Image
    401             try:
    402                 from PIL import PILLOW_VERSION
    403                 self.main.PackageVersions.append([Image.__name__,PILLOW_VERSION])
    404             except:
    405                 self.main.PackageVersions.append([Image.__name__,Image.VERSION])
    406         except ImportError:
    407             try:
    408                 import Image
    409                 self.main.PackageVersions.append([Image.__name__,Image.VERSION])
    410             except ImportError:
    411                 pass
    412         self.main.PackageVersions.append(['Platform',sys.platform+' '+platform.architecture()[0]+' '+platform.machine()])
     395        self.main.PackageVersions = G2fil.get_python_versions([wx, mpl, np, sp, ogl])
    413396       
    414397        return True
    … …  
    578561        '''
    579562
    580         self.ImportPhaseReaderlist = []
    581         self._init_Import_routines('phase',self.ImportPhaseReaderlist,'Phase')
    582         self.ImportSfactReaderlist = []
    583         self._init_Import_routines('sfact',self.ImportSfactReaderlist,'Struct_Factor')
    584         self.ImportPowderReaderlist = []
    585         self._init_Import_routines('pwd',self.ImportPowderReaderlist,'Powder_Data')
    586         self.ImportSmallAngleReaderlist = []
    587         self._init_Import_routines('sad',self.ImportSmallAngleReaderlist,'SmallAngle_Data')
    588         self.ImportReflectometryReaderlist = []
    589         self._init_Import_routines('rfd',self.ImportReflectometryReaderlist,'Reflectometry_Data')
    590         self.ImportPDFReaderlist = []
    591         self._init_Import_routines('pdf',self.ImportPDFReaderlist,'PDF_Data')
    592         self.ImportImageReaderlist = []
    593         self._init_Import_routines('img',self.ImportImageReaderlist,'Images')
     563        self.ImportPhaseReaderlist = G2fil.LoadImportRoutines('phase','Phase')
     564        self.ImportSfactReaderlist = G2fil.LoadImportRoutines('sfact','Struct_Factor')
     565        self.ImportPowderReaderlist = G2fil.LoadImportRoutines('pwd','Powder_Data')
     566        self.ImportSmallAngleReaderlist = G2fil.LoadImportRoutines('sad','SmallAngle_Data')
     567        self.ImportReflectometryReaderlist = G2fil.LoadImportRoutines('rfd','Reflectometry_Data')
     568        self.ImportPDFReaderlist = G2fil.LoadImportRoutines('pdf','PDF_Data')
     569        self.ImportImageReaderlist = G2fil.LoadImportRoutines('img','Images')
    594570        self.ImportMenuId = {}
    595 
    596     def _init_Import_routines(self,prefix,readerlist,errprefix):
    597         '''import all the import readers matching the prefix
    598         '''
    599         #path2GSAS2 = os.path.dirname(os.path.realpath(__file__)) # location of this file
    600         #pathlist = sys.path[:]
    601         #if path2GSAS2 not in pathlist: pathlist.append(path2GSAS2)
    602         #path2GSAS2 = os.path.join(
    603         #    os.path.dirname(os.path.realpath(__file__)), # location of this file
    604         #    'imports')
    605         pathlist = sys.path[:]
    606         #if path2GSAS2 not in pathlist: pathlist.append(path2GSAS2)
    607         if '.' not in pathlist: pathlist.append('.') # insert the directory where G2 is started
    608 
    609         filelist = []
    610         for path in pathlist:
    611             for filename in glob.iglob(os.path.join(
    612                 path,
    613                 "G2"+prefix+"*.py")):
    614                 filelist.append(filename)   
    615                 #print 'debug: found',filename
    616         filelist = sorted(list(set(filelist))) # remove duplicates
    617         for filename in filelist:
    618             path,rootname = os.path.split(filename)
    619             pkg = os.path.splitext(rootname)[0]
    620             try:
    621                 fp = None
    622                 fp, fppath,desc = imp.find_module(pkg,[path,])
    623                 pkg = imp.load_module(pkg,fp,fppath,desc)
    624                 for clss in inspect.getmembers(pkg): # find classes defined in package
    625                     if clss[0].startswith('_'): continue
    626                     if inspect.isclass(clss[1]):
    627                         # check if we have the required methods
    628                         for m in 'Reader','ExtensionValidator','ContentsValidator':
    629                             if not hasattr(clss[1],m): break
    630                             if not callable(getattr(clss[1],m)): break
    631                         else:
    632                             reader = clss[1]() # create an import instance
    633                             if reader.UseReader:
    634                                 readerlist.append(reader)
    635             except AttributeError:
    636                 print 'Import_'+errprefix+': Attribute Error '+ filename
    637             #except ImportError:
    638             #    print 'Import_'+errprefix+': Error importing file '+ filename
    639             except Exception,errmsg:
    640                 print('\nImport_'+errprefix+': Error importing file '+ filename)
    641                 print(u'Error message: {}\n'.format(errmsg))
    642             if fp: fp.close()
    643571
    644572    def EnableSeqRefineMenu(self):
    … …  
    13111239        finds matching bank no. to load - problem with nonmatches?
    13121240
     1241        Note that this routine performs a similar role to :func:`GSASIIfiles.ReadPowderInstprm`,
     1242        but this will call a GUI routine for selection when needed. TODO: refactor to combine
     1243
    13131244        :param list instLines: strings from GSAS-II parameter file; can be concatenated with ';'
    13141245        :param int  bank: bank number to check when instprm file has '#BANK n:...' strings
    … …  
    13911322            newVals.append(val)
    13921323            il += 1
    1393         return [G2IO.makeInstDict(newItems,newVals,len(newVals)*[False,]),{}]
     1324        return [G2fil.makeInstDict(newItems,newVals,len(newVals)*[False,]),{}]
    13941325       
    13951326    def ReadPowderIparm(self,instfile,bank,databanks,rd):
    … …  
    14661397        :returns: a list of two dicts, the first containing instrument parameters
    14671398          and the second used for TOF lookup tables for profile coeff.
    1468 
    14691399        '''
    1470         def SetPowderInstParms(Iparm, rd):
    1471             '''extracts values from instrument parameters in rd.instdict
    1472             or in array Iparm.
    1473             Create and return the contents of the instrument parameter tree entry.
    1474             '''
    1475             Irads = {0:' ',1:'CrKa',2:'FeKa',3:'CuKa',4:'MoKa',5:'AgKa',6:'TiKa',7:'CoKa'}
    1476             DataType = Iparm['INS   HTYPE '].strip()[:3]  # take 1st 3 chars
    1477             # override inst values with values read from data file
    1478             Bank = rd.powderentry[2]    #should be used in multibank iparm files
    1479             if rd.instdict.get('type'):
    1480                 DataType = rd.instdict.get('type')
    1481             data = [DataType,]
    1482             instname = Iparm.get('INS  1INAME ')
    1483             irad = int(Iparm.get('INS  1 IRAD ','0'))
    1484             if instname:
    1485                 rd.Sample['InstrName'] = instname.strip()
    1486             if 'C' in DataType:
    1487                 wave1 = None
    1488                 wave2 = 0.0
    1489                 if rd.instdict.get('wave'):
    1490                     wl = rd.instdict.get('wave')
    1491                     wave1 = wl[0]
    1492                     if len(wl) > 1: wave2 = wl[1]
    1493                 s = Iparm['INS  1 ICONS']
    1494                 if not wave1:
    1495                     wave1 = G2IO.sfloat(s[:10])
    1496                     wave2 = G2IO.sfloat(s[10:20])
    1497                 v = (wave1,wave2,
    1498                      G2IO.sfloat(s[20:30])/100.,G2IO.sfloat(s[55:65]),G2IO.sfloat(s[40:50])) #get lam1, lam2, zero, pola & ratio
    1499                 if not v[1]:
    1500                     names = ['Type','Lam','Zero','Polariz.','U','V','W','X','Y','SH/L','Azimuth']
    1501                     v = (v[0],v[2],v[4])
    1502                     codes = [0,0,0,0]
    1503                     rd.Sample.update({'Type':'Debye-Scherrer','Absorption':[0.,False],'DisplaceX':[0.,False],'DisplaceY':[0.,False]})
    1504                 else:
    1505                     names = ['Type','Lam1','Lam2','Zero','I(L2)/I(L1)','Polariz.','U','V','W','X','Y','SH/L','Azimuth']
    1506                     codes = [0,0,0,0,0,0]
    1507                     rd.Sample.update({'Type':'Bragg-Brentano','Shift':[0.,False],'Transparency':[0.,False],
    1508                         'SurfRoughA':[0.,False],'SurfRoughB':[0.,False]})
    1509                 data.extend(v)
    1510                 if 'INS  1PRCF  ' in Iparm:
    1511                     v1 = Iparm['INS  1PRCF  '].split()                                                 
    1512                     v = Iparm['INS  1PRCF 1'].split()
    1513                     data.extend([float(v[0]),float(v[1]),float(v[2])])                  #get GU, GV & GW - always here
    1514                     azm = float(Iparm.get('INS  1DETAZM','0.0'))
    1515                     v = Iparm['INS  1PRCF 2'].split()
    1516                     if v1[0] == 3:
    1517                         data.extend([float(v[0]),float(v[1]),float(v[2])+float(v[3],azm)])  #get LX, LY, S+H/L & azimuth
    1518                     else:
    1519                         data.extend([0.0,0.0,0.002,azm])                                      #OK defaults if fxn #3 not 1st in iprm file                   
    1520                 else:
    1521                     v1 = Iparm['INS  1PRCF1 '].split()                                                 
    1522                     v = Iparm['INS  1PRCF11'].split()
    1523                     data.extend([float(v[0]),float(v[1]),float(v[2])])                  #get GU, GV & GW - always here
    1524                     azm = float(Iparm.get('INS  1DETAZM','0.0'))
    1525                     v = Iparm['INS  1PRCF12'].split()
    1526                     if v1[0] == 3:
    1527                         data.extend([float(v[0]),float(v[1]),float(v[2])+float(v[3],azm)])  #get LX, LY, S+H/L & azimuth
    1528                     else:
    1529                         data.extend([0.0,0.0,0.002,azm])                                      #OK defaults if fxn #3 not 1st in iprm file
    1530                 codes.extend([0,0,0,0,0,0,0])
    1531                 Iparm1 = G2IO.makeInstDict(names,data,codes)
    1532                 Iparm1['Source'] = [Irads[irad],Irads[irad]]
    1533                 Iparm1['Bank'] = [Bank,Bank,0]
    1534                 return [Iparm1,{}]
    1535             elif 'T' in DataType:
    1536                 names = ['Type','fltPath','2-theta','difC','difA', 'difB','Zero','alpha','beta-0','beta-1',
    1537                     'beta-q','sig-0','sig-1','sig-2','sig-q', 'X','Y','Azimuth',]
    1538                 codes = [0,0,0,0,0, 0,0,0,0,0, 0,0,0,0,0, 0,0,0,]
    1539                 azm = 0.
    1540                 if 'INS  1DETAZM' in Iparm:
    1541                     azm = float(Iparm['INS  1DETAZM'])
    1542                 rd.Sample['Azimuth'] = azm
    1543                 fltPath0 = 20.                      #arbitrary
    1544                 if 'INS   FPATH1' in Iparm:                   
    1545                     s = Iparm['INS   FPATH1'].split()
    1546                     fltPath0 = G2IO.sfloat(s[0])
    1547                 if 'INS  1BNKPAR' not in Iparm:     #bank missing from Iparm file
    1548                     return []
    1549                 s = Iparm['INS  1BNKPAR'].split()
    1550                 fltPath1 = G2IO.sfloat(s[0])
    1551                 data.extend([fltPath0+fltPath1,])               #Flight path source-sample-detector
    1552                 data.extend([G2IO.sfloat(s[1]),])               #2-theta for bank
    1553                 s = Iparm['INS  1 ICONS'].split()
    1554                 data.extend([G2IO.sfloat(s[0]),G2IO.sfloat(s[1]),0.0,G2IO.sfloat(s[2])])    #difC,difA,difB,Zero
    1555                 if 'INS  1PRCF  ' in Iparm:
    1556                     s = Iparm['INS  1PRCF  '].split()
    1557                     pfType = int(s[0])
    1558                     s = Iparm['INS  1PRCF 1'].split()
    1559                     if abs(pfType) == 1:
    1560                         data.extend([G2IO.sfloat(s[1]),G2IO.sfloat(s[2]),G2IO.sfloat(s[3])]) #alpha, beta-0, beta-1
    1561                         s = Iparm['INS  1PRCF 2'].split()
    1562                         data.extend([0.0,0.0,G2IO.sfloat(s[1]),G2IO.sfloat(s[2]),0.0,0.0,0.0,azm])    #beta-q, sig-0, sig-1, sig-2, sig-q, X, Y
    1563                     elif abs(pfType) in [3,4,5]:
    1564                         data.extend([G2IO.sfloat(s[0]),G2IO.sfloat(s[1]),G2IO.sfloat(s[2])]) #alpha, beta-0, beta-1
    1565                         if abs(pfType) == 4:
    1566                             data.extend([0.0,0.0,G2IO.sfloat(s[3]),0.0,0.0,0.0,0.0,azm])    #beta-q, sig-0, sig-1, sig-2, sig-q, X, Y
    1567                         else:
    1568                             s = Iparm['INS  1PRCF 2'].split()
    1569                             data.extend([0.0,0.0,G2IO.sfloat(s[0]),G2IO.sfloat(s[1]),0.0,0.0,0.0,azm])    #beta-q, sig-0, sig-1, sig-2, sig-q, X, Y                       
    1570                     elif abs(pfType) == 2:
    1571                         data.extend([G2IO.sfloat(s[1]),0.0,1./G2IO.sfloat(s[3])]) #alpha, beta-0, beta-1
    1572                         data.extend([0.0,0.0,G2IO.sfloat(s[1]),0.0,0.0,0.0,0.0,azm])    #beta-q, sig-0, sig-1, sig-2, sig-q, X, Y                           
    1573                 else:
    1574                     s = Iparm['INS  1PRCF1 '].split()
    1575                     pfType = int(s[0])
    1576                     s = Iparm['INS  1PRCF11'].split()
    1577                     if abs(pfType) == 1:
    1578                         data.extend([G2IO.sfloat(s[1]),G2IO.sfloat(s[2]),G2IO.sfloat(s[3])]) #alpha, beta-0, beta-1
    1579                         s = Iparm['INS  1PRCF12'].split()
    1580                         data.extend([0.0,0.0,G2IO.sfloat(s[1]),G2IO.sfloat(s[2]),0.0,0.0,0.0,azm])    #beta-q, sig-0, sig-1, sig-2, sig-q, X, Y
    1581                     elif abs(pfType) in [3,4,5]:
    1582                         data.extend([G2IO.sfloat(s[0]),G2IO.sfloat(s[1]),G2IO.sfloat(s[2])]) #alpha, beta-0, beta-1
    1583                         if abs(pfType) == 4:
    1584                             data.extend([0.0,0.0,G2IO.sfloat(s[3]),0.0,0.0,0.0,0.0,azm])    #beta-q, sig-0, sig-1, sig-2, sig-q, X, Y
    1585                         else:
    1586                             s = Iparm['INS  1PRCF12'].split()
    1587                             data.extend([0.0,0.0,G2IO.sfloat(s[0]),G2IO.sfloat(s[1]),0.0,0.0,0.0,azm])    #beta-q, sig-0, sig-1, sig-2, sig-q, X, Y                       
    1588                 Inst1 = G2IO.makeInstDict(names,data,codes)
    1589                 Inst1['Bank'] = [Bank,Bank,0]
    1590                 Inst2 = {}
    1591                 if pfType < 0:
    1592                     Ipab = 'INS  1PAB'+str(-pfType)
    1593                     Npab = int(Iparm[Ipab+'  '].strip())
    1594                     Inst2['Pdabc'] = []
    1595                     for i in range(Npab):
    1596                         k = Ipab+str(i+1).rjust(2)
    1597                         s = Iparm[k].split()
    1598                         Inst2['Pdabc'].append([float(t) for t in s])
    1599                     Inst2['Pdabc'] = np.array(Inst2['Pdabc'])
    1600                     Inst2['Pdabc'].T[3] += Inst2['Pdabc'].T[0]*Inst1['difC'][0] #turn 3rd col into TOF
    1601                 if 'INS  1I ITYP' in Iparm:
    1602                     s = Iparm['INS  1I ITYP'].split()
    1603                     Ityp = int(s[0])
    1604                     Tminmax = [float(s[1])*1000.,float(s[2])*1000.]
    1605                     Itypes = ['Exponential','Maxwell/Exponential','','Maxwell/Chebyschev','']
    1606                     if Ityp in [1,2,4]:
    1607                         Inst2['Itype'] = Itypes[Ityp-1]
    1608                         Inst2['Tminmax'] = Tminmax
    1609                         Icoeff = []
    1610                         Iesd = []
    1611                         Icovar = []                   
    1612                         for i in range(3):
    1613                             s = Iparm['INS  1ICOFF'+str(i+1)].split()
    1614                             Icoeff += [float(S) for S in s]
    1615                             s = Iparm['INS  1IECOF'+str(i+1)].split()
    1616                             Iesd += [float(S) for S in s]
    1617                         NT = 10
    1618                         for i in range(8):
    1619                             s = Iparm['INS  1IECOR'+str(i+1)]
    1620                             if i == 7:
    1621                                 NT = 8
    1622                             Icovar += [float(s[6*j:6*j+6]) for j in range(NT)]
    1623                         Inst2['Icoeff'] = Icoeff
    1624                         Inst2['Iesd'] = Iesd
    1625                         Inst2['Icovar'] = Icovar
    1626                 return [Inst1,Inst2]
    16271400               
    16281401        def GetDefaultParms(self,rd):
    … …  
    17041477                    rd.instfile = instfile
    17051478                    rd.instmsg = instfile + ' bank ' + str(rd.instbank)
    1706                     return SetPowderInstParms(Iparm,rd)
     1479                    return G2fil.SetPowderInstParms(Iparm,rd)
    17071480                else:
    17081481                    #print 'debug: open/read failed',instfile
    … …  
    17331506                        rd.instfile = instfile
    17341507                        rd.instmsg = instfile + ' bank ' + str(rd.instbank)
    1735                         instParmList = SetPowderInstParms(Iparm,rd)     #this is [Inst1,Inst2] a pair of dicts
     1508                        instParmList = G2fil.SetPowderInstParms(Iparm,rd)     #this is [Inst1,Inst2] a pair of dicts
    17361509                if 'list' in str(type(instParmList)):   #record stuff & return stuff
    17371510                    rd.instfile = instfile
    … …  
    17791552                    rd.instfile = instfile
    17801553                    rd.instmsg = instfile + ' bank ' + str(rd.instbank)
    1781                     return SetPowderInstParms(Iparm,rd)
     1554                    return G2fil.SetPowderInstParms(Iparm,rd)
    17821555                else:
    17831556                    self.ErrorDialog('Read Error',
    … …  
    31872960                names = ['Type','Lam','Zero']
    31882961                codes = [0,0,0]
    3189                 inst = [G2IO.makeInstDict(names,data,codes),{}]
     2962                inst = [G2fil.makeInstDict(names,data,codes),{}]
    31902963                self.GPXtree.SetItemPyData(self.GPXtree.AppendItem(Id,text='Instrument Parameters'),inst)
    31912964                self.GPXtree.SetItemPyData(self.GPXtree.AppendItem(Id,text='Comments'),comments)
    … …  
    73577130    G2frame.dataWindow.currentGrids = []
    73587131    G2frame.dataDisplay = G2G.GSGrid(parent=G2frame.dataWindow)
     7132    G2frame.dataDisplay.SetScrollRate(1,1)
    73597133    G2frame.dataWindow.GetSizer().Add(G2frame.dataDisplay,1,wx.ALL|wx.EXPAND)
    73607134    G2frame.SeqTable = G2G.Table([list(cl) for cl in zip(*G2frame.colList)],     # convert from columns to rows

  • branch/2frame/GSASIIpwdGUI.py

    r2960 r2963  
    3232import GSASIIpwd as G2pwd
    3333import GSASIIIO as G2IO
     34import GSASIIfiles as G2fil
    3435import GSASIIobj as G2obj
    3536import GSASIIlattice as G2lat
    … …  
    14701471                    if 'Bank' not in Inst:  #patch for old .instprm files - may cause faults for TOF data
    14711472                        Inst['Bank'] = [1,1,0]
    1472                     data = G2IO.makeInstDict(newItems,newVals,len(newVals)*[False,])
     1473                    data = G2fil.makeInstDict(newItems,newVals,len(newVals)*[False,])
    14731474                    G2frame.GPXtree.SetItemPyData(G2gd.GetGPXtreeItemId(G2frame,G2frame.PatternId,'Instrument Parameters'),[data,Inst2])
    14741475                    RefreshInstrumentGrid(event,doAnyway=True)          #to get peaks updated

  • branch/2frame/GSASIIscriptable.py

    r2895 r2963  
    1010
    1111"""
     12from __future__ import division, print_function # needed?
    1213import os.path as ospath
    1314import sys
    … …  
    1516import imp
    1617import copy
     18import platform
     19import os
     20import random as ran
     21import json
     22
     23import numpy.ma as ma
     24import scipy.interpolate as si
     25import numpy as np
     26import scipy as sp
     27import matplotlib as mpl
     28
    1729import GSASIIpath
     30GSASIIpath.SetBinaryPath() # would rather have this in __name__ == '__main__' stanza
     31import GSASIIIO as G2IO
     32import GSASIIfiles as G2fil
    1833import GSASIIobj as G2obj
     34import GSASIIpwd as G2pwd
     35import GSASIIstrIO as G2strIO
     36import GSASIIspc as G2spc
     37import GSASIIstrMain as G2strMain
     38import GSASIImath as G2mth
     39import GSASIIElem as G2elem
    1940
    2041def LoadDictFromProjFile(ProjFile):
    21     ''''Read a GSAS-II project file and load items to dictionary
     42    '''Read a GSAS-II project file and load items to dictionary
    2243    :param: str ProjFile: GSAS-II project (name.gpx) full file name
    2344    :returns: dict Project: representation of gpx file following the GSAS-II tree struture
    … …  
    2546        data dict = {'data':item data whch may be list, dict or None,'subitems':subdata (if any)}
    2647    :returns: list nameList: names of main tree entries & subentries used to reconstruct project file
     48
    2749    Example for fap.gpx:
    2850    Project = {                 #NB:dict order is not tree order
    … …  
    6082        return
    6183    file = open(ProjFile,'rb')
    62     print 'loading from file: ',ProjFile
     84    print('loading from file: {}'.format(ProjFile))
    6385    Project = {}
    6486    nameList = []
    … …  
    90112    '''
    91113    file = open(ProjFile,'wb')
    92     print 'save to file: ',ProjFile
     114    print('save to file: {}'.format(ProjFile))
    93115    for name in nameList:
    94116        data = []
    … …  
    118140    rd = rdclass.GSAS_ReaderClass()   
    119141    if not rd.scriptable:
    120         print '**** ERROR: ',reader,' is not a scriptable reader'
     142        print(u'**** ERROR: '+reader+u' is not a scriptable reader')
    121143        return None
    122144    fl = open(filename,'rb')
    … …  
    137159                    repeat = True
    138160        return rdlist
    139     print rd.errors
     161    print(rd.errors)
    140162    return None
    141163   
     164def SetDefaultDData(dType,histoName,NShkl=0,NDij=0):
     165    '''Create an initial Histogram dictionary
     166   
     167    Author: Jackson O'Donnell jacksonhodonnell@gmail.com
     168    '''
     169    if dType in ['SXC','SNC']:
     170        return {'Histogram':histoName,'Show':False,'Scale':[1.0,True],
     171            'Babinet':{'BabA':[0.0,False],'BabU':[0.0,False]},
     172            'Extinction':['Lorentzian','None', {'Tbar':0.1,'Cos2TM':0.955,
     173            'Eg':[1.e-10,False],'Es':[1.e-10,False],'Ep':[1.e-10,False]}],
     174            'Flack':[0.0,False]}
     175    elif dType == 'SNT':
     176        return {'Histogram':histoName,'Show':False,'Scale':[1.0,True],
     177            'Babinet':{'BabA':[0.0,False],'BabU':[0.0,False]},
     178            'Extinction':['Lorentzian','None', {
     179            'Eg':[1.e-10,False],'Es':[1.e-10,False],'Ep':[1.e-10,False]}]}
     180    elif 'P' in dType:
     181        return {'Histogram':histoName,'Show':False,'Scale':[1.0,False],
     182            'Pref.Ori.':['MD',1.0,False,[0,0,1],0,{},[],0.1],
     183            'Size':['isotropic',[1.,1.,1.],[False,False,False],[0,0,1],
     184                [1.,1.,1.,0.,0.,0.],6*[False,]],
     185            'Mustrain':['isotropic',[1000.0,1000.0,1.0],[False,False,False],[0,0,1],
     186                NShkl*[0.01,],NShkl*[False,]],
     187            'HStrain':[NDij*[0.0,],NDij*[False,]],
     188            'Extinction':[0.0,False],'Babinet':{'BabA':[0.0,False],'BabU':[0.0,False]}}
     189
     190
     191def PreSetup(data):
     192    '''Create part of an initial (empty) phase dictionary
     193   
     194    from GSASIIphsGUI.py, near end of UpdatePhaseData
     195    Author: Jackson O'Donnell jacksonhodonnell@gmail.com
     196    '''
     197    if 'RBModels' not in data:
     198        data['RBModels'] = {}
     199    if 'MCSA' not in data:
     200        data['MCSA'] = {'Models':[{'Type':'MD','Coef':[1.0,False,[.8,1.2],],'axis':[0,0,1]}],'Results':[],'AtInfo':{}}
     201    if 'dict' in str(type(data['MCSA']['Results'])):
     202        data['MCSA']['Results'] = []
     203    if 'Modulated' not in data['General']:
     204        data['General']['Modulated'] = False
     205    if 'modulated' in data['General']['Type']:
     206        data['General']['Modulated'] = True
     207        data['General']['Type'] = 'nuclear'
     208
     209
     210def SetupGeneral(data, dirname):
     211    """Helps initialize phase data.
     212
     213    From GSASIIphsGui.py, function of the same name. Minor changes for imports etc.
     214    Author: Jackson O'Donnell jacksonhodonnell@gmail.com
     215    """
     216    mapDefault = {'MapType':'','RefList':'','Resolution':0.5,'Show bonds':True,
     217                'rho':[],'rhoMax':0.,'mapSize':10.0,'cutOff':50.,'Flip':False}
     218    generalData = data['General']
     219    atomData = data['Atoms']
     220    generalData['AtomTypes'] = []
     221    generalData['Isotopes'] = {}
     222
     223    if 'Isotope' not in generalData:
     224        generalData['Isotope'] = {}
     225    if 'Data plot type' not in generalData:
     226        generalData['Data plot type'] = 'Mustrain'
     227    if 'POhkl' not in generalData:
     228        generalData['POhkl'] = [0,0,1]
     229    if 'Map' not in generalData:
     230        generalData['Map'] = mapDefault.copy()
     231    if 'Flip' not in generalData:
     232        generalData['Flip'] = {'RefList':'','Resolution':0.5,'Norm element':'None',
     233            'k-factor':0.1,'k-Max':20.,}
     234    if 'testHKL' not in generalData['Flip']:
     235        generalData['Flip']['testHKL'] = [[0,0,2],[2,0,0],[1,1,1],[0,2,0],[1,2,3]]
     236    if 'doPawley' not in generalData:
     237        generalData['doPawley'] = False     #ToDo: change to ''
     238    if 'Pawley dmin' not in generalData:
     239        generalData['Pawley dmin'] = 1.0
     240    if 'Pawley neg wt' not in generalData:
     241        generalData['Pawley neg wt'] = 0.0
     242    if 'Algolrithm' in generalData.get('MCSA controls',{}) or \
     243        'MCSA controls' not in generalData:
     244        generalData['MCSA controls'] = {'Data source':'','Annealing':[50.,0.001,50],
     245        'dmin':2.0,'Algorithm':'log','Jump coeff':[0.95,0.5],'boltzmann':1.0,
     246        'fast parms':[1.0,1.0,1.0],'log slope':0.9,'Cycles':1,'Results':[],'newDmin':True}
     247    if 'AtomPtrs' not in generalData:
     248        generalData['AtomPtrs'] = [3,1,7,9]
     249        if generalData['Type'] == 'macromolecular':
     250            generalData['AtomPtrs'] = [6,4,10,12]
     251        elif generalData['Type'] == 'magnetic':
     252            generalData['AtomPtrs'] = [3,1,10,12]
     253    if generalData['Type'] in ['modulated',]:
     254        generalData['Modulated'] = True
     255        generalData['Type'] = 'nuclear'
     256        if 'Super' not in generalData:
     257            generalData['Super'] = 1
     258            generalData['SuperVec'] = [[0,0,.1],False,4]
     259            generalData['SSGData'] = {}
     260        if '4DmapData' not in generalData:
     261            generalData['4DmapData'] = mapDefault.copy()
     262            generalData['4DmapData'].update({'MapType':'Fobs'})
     263    if 'Modulated' not in generalData:
     264        generalData['Modulated'] = False
     265    if 'HydIds' not in generalData:
     266        generalData['HydIds'] = {}
     267    cx,ct,cs,cia = generalData['AtomPtrs']
     268    generalData['NoAtoms'] = {}
     269    generalData['BondRadii'] = []
     270    generalData['AngleRadii'] = []
     271    generalData['vdWRadii'] = []
     272    generalData['AtomMass'] = []
     273    generalData['Color'] = []
     274    if generalData['Type'] == 'magnetic':
     275        generalData['MagDmin'] = generalData.get('MagDmin',1.0)
     276        landeg = generalData.get('Lande g',[])
     277    generalData['Mydir'] = dirname
     278    badList = {}
     279    for iat,atom in enumerate(atomData):
     280        atom[ct] = atom[ct].lower().capitalize()              #force to standard form
     281        if generalData['AtomTypes'].count(atom[ct]):
     282            generalData['NoAtoms'][atom[ct]] += atom[cx+3]*float(atom[cs+1])
     283        elif atom[ct] != 'UNK':
     284            Info = G2elem.GetAtomInfo(atom[ct])
     285            if not Info:
     286                if atom[ct] not in badList:
     287                    badList[atom[ct]] = 0
     288                badList[atom[ct]] += 1
     289                atom[ct] = 'UNK'
     290                continue
     291            atom[ct] = Info['Symbol'] # N.B. symbol might be changed by GetAtomInfo
     292            generalData['AtomTypes'].append(atom[ct])
     293            generalData['Z'] = Info['Z']
     294            generalData['Isotopes'][atom[ct]] = Info['Isotopes']
     295            generalData['BondRadii'].append(Info['Drad'])
     296            generalData['AngleRadii'].append(Info['Arad'])
     297            generalData['vdWRadii'].append(Info['Vdrad'])
     298            if atom[ct] in generalData['Isotope']:
     299                if generalData['Isotope'][atom[ct]] not in generalData['Isotopes'][atom[ct]]:
     300                    isotope = generalData['Isotopes'][atom[ct]].keys()[-1]
     301                    generalData['Isotope'][atom[ct]] = isotope
     302                generalData['AtomMass'].append(Info['Isotopes'][generalData['Isotope'][atom[ct]]]['Mass'])
     303            else:
     304                generalData['Isotope'][atom[ct]] = 'Nat. Abund.'
     305                if 'Nat. Abund.' not in generalData['Isotopes'][atom[ct]]:
     306                    isotope = generalData['Isotopes'][atom[ct]].keys()[-1]
     307                    generalData['Isotope'][atom[ct]] = isotope
     308                generalData['AtomMass'].append(Info['Mass'])
     309            generalData['NoAtoms'][atom[ct]] = atom[cx+3]*float(atom[cs+1])
     310            generalData['Color'].append(Info['Color'])
     311            if generalData['Type'] == 'magnetic':
     312                if len(landeg) < len(generalData['AtomTypes']):
     313                    landeg.append(2.0)
     314    if generalData['Type'] == 'magnetic':
     315        generalData['Lande g'] = landeg[:len(generalData['AtomTypes'])]
     316
     317    if badList:
     318        msg = 'Warning: element symbol(s) not found:'
     319        for key in badList:
     320            msg += '\n\t' + key
     321            if badList[key] > 1:
     322                msg += ' (' + str(badList[key]) + ' times)'
     323        raise Exception("Phase error:\n" + msg)
     324        # wx.MessageBox(msg,caption='Element symbol error')
     325    F000X = 0.
     326    F000N = 0.
     327    for i,elem in enumerate(generalData['AtomTypes']):
     328        F000X += generalData['NoAtoms'][elem]*generalData['Z']
     329        isotope = generalData['Isotope'][elem]
     330        F000N += generalData['NoAtoms'][elem]*generalData['Isotopes'][elem][isotope]['SL'][0]
     331    generalData['F000X'] = F000X
     332    generalData['F000N'] = F000N
     333    generalData['Mass'] = G2mth.getMass(generalData)
     334
     335def make_empty_project(author=None, filename=None):
     336    """Creates an dictionary in the style of GSASIIscriptable, for an empty
     337    project.
     338
     339    If no author name or filename is supplied, 'no name' and
     340    <current dir>/test_output.gpx are used , respectively.
     341
     342    Returns: project dictionary, name list
     343    Author: Jackson O'Donnell jacksonhodonnell@gmail.com
     344    """
     345    if not filename:
     346        filename = os.path.join(os.getcwd(), 'test_output.gpx')
     347    else:
     348        filename = os.path.abspath(filename)
     349    gsasii_version = str(GSASIIpath.GetVersionNumber())
     350    python_library_versions = G2fil.get_python_versions([mpl, np, sp])
     351
     352    controls_data = dict(G2obj.DefaultControls)
     353    controls_data['LastSavedAs'] = unicode(filename)
     354    controls_data['LastSavedUsing'] = gsasii_version
     355    controls_data['PythonVersions'] = python_library_versions
     356    if author:
     357        controls_data['Author'] = author
     358
     359    output = {'Constraints': {'data': {'HAP': [], 'Hist': [], 'Phase': [],
     360                                       'Global': []}},
     361              'Controls': {'data': controls_data},
     362              u'Covariance': {'data': {}},
     363              u'Notebook': {'data': ['']},
     364              u'Restraints': {'data': {}},
     365              u'Rigid bodies': {'data': {'RBIds': {'Residue': [], 'Vector': []},
     366                                'Residue': {'AtInfo': {}},
     367                                'Vector':  {'AtInfo': {}}}}}
     368
     369    names = [[u'Notebook'], [u'Controls'], [u'Covariance'],
     370             [u'Constraints'], [u'Restraints'], [u'Rigid bodies']]
     371
     372    return output, names
     373
     374
     375class G2ImportException(Exception):
     376    pass
     377
     378
     379def import_generic(filename, readerlist):
     380    """Attempt to import a filename, using a list of reader objects.
     381
     382    Returns the first reader object which worked."""
     383    # Translated from OnImportGeneric method in GSASII.py
     384    primaryReaders, secondaryReaders = [], []
     385    for reader in readerlist:
     386        flag = reader.ExtensionValidator(filename)
     387        if flag is None:
     388            secondaryReaders.append(reader)
     389        elif flag:
     390            primaryReaders.append(reader)
     391    if not secondaryReaders and not primaryReaders:
     392        raise G2ImportException("Could not read file: ", filename)
     393
     394    with open(filename, 'Ur') as fp:
     395        rd_list = []
     396
     397        for rd in primaryReaders + secondaryReaders:
     398            # Initialize reader
     399            rd.selections = []
     400            rd.dnames = []
     401            rd.ReInitialize()
     402            # Rewind file
     403            fp.seek(0)
     404            rd.errors = ""
     405            if not rd.ContentsValidator(fp):
     406                # Report error
     407                pass
     408            if len(rd.selections) > 1:
     409                # Select data?
     410                # GSASII.py:543
     411                raise G2ImportException("Not sure what data to select")
     412
     413            block = 0
     414            rdbuffer = {}
     415            fp.seek(0)
     416            repeat = True
     417            while repeat:
     418                repeat = False
     419                block += 1
     420                rd.objname = os.path.basename(filename)
     421                flag = rd.Reader(filename, fp, buffer=rdbuffer, blocknum=block)
     422                if flag:
     423                    # Omitting image loading special cases
     424                    rd.readfilename = filename
     425                    rd_list.append(copy.deepcopy(rd))
     426                    repeat = rd.repeat
     427                else:
     428                    raise G2ImportException("{}.Reader() returned:".format(rd),
     429                                            flag)
     430
     431            if rd_list:
     432                if rd.warnings:
     433                    print("Read warning by", rd.formatName, "reader:",
     434                          rd.warnings, file=sys.stderr)
     435                return rd_list
     436    raise G2ImportException("No reader could read file: " + filename)
     437
     438
     439def load_iprms(instfile, reader):
     440    """Loads instrument parameters from a file, and edits the
     441    given reader.
     442
     443    Returns a 2-tuple of (Iparm1, Iparm2) parameters
     444    """
     445    ext = os.path.splitext(instfile)[1]
     446
     447    if ext.lower() == '.instprm':
     448        # New GSAS File, load appropriately
     449        with open(instfile) as f:
     450            lines = f.readlines()
     451        bank = reader.powderentry[2]
     452        numbanks = reader.numbanks
     453        iparms = G2fil.ReadPowderInstprm(lines, bank, numbanks, reader)
     454
     455        reader.instfile = instfile
     456        reader.instmsg = 'GSAS-II file' + instfile
     457        return iparms
     458    elif ext.lower() not in ('.prm', '.inst', '.ins'):
     459        raise ValueError('Expected .prm file, found: ', instfile)
     460
     461    # It's an old GSAS file, load appropriately
     462    Iparm = {}
     463    with open(instfile, 'Ur') as fp:
     464        for line in fp:
     465            if '#' in line:
     466                continue
     467            Iparm[line[:12]] = line[12:-1]
     468    ibanks = int(Iparm.get('INS   BANK  ', '1').strip())
     469    if ibanks == 1:
     470        reader.instbank = 1
     471        reader.powderentry[2] = 1
     472        reader.instfile = instfile
     473        reader.instmsg = instfile + ' bank ' + str(reader.instbank)
     474        return G2fil.SetPowderInstParms(Iparm, reader)
     475    # TODO handle >1 banks
     476    raise NotImplementedError("Check GSASIIfiles.py: ReadPowderInstprm")
     477
     478def load_pwd_from_reader(reader, instprm, existingnames=[]):
     479    """Loads powder data from a reader object, and assembles it into a GSASII data tree.
     480
     481    Returns: (name, tree) - 2-tuple of the histogram name (str), and data
     482    Author: Jackson O'Donnell jacksonhodonnell@gmail.com
     483    """
     484    HistName = 'PWDR ' + G2obj.StripUnicode(reader.idstring, '_')
     485    HistName = unicode(G2obj.MakeUniqueLabel(HistName, existingnames))
     486
     487    try:
     488        Iparm1, Iparm2 = instprm
     489    except ValueError:
     490        # TODO load_iprms
     491        Iparm1, Iparm2 = load_iprms(instprm, reader)
     492
     493    Ymin = np.min(reader.powderdata[1])
     494    Ymax = np.max(reader.powderdata[1])
     495    valuesdict = {'wtFactor': 1.0,
     496                  'Dummy': False,
     497                  'ranId': ran.randint(0, sys.maxint),
     498                  'Offset': [0.0, 0.0], 'delOffset': 0.02*Ymax,
     499                  'refOffset': -0.1*Ymax, 'refDelt': 0.1*Ymax,
     500                  'qPlot': False, 'dPlot': False, 'sqrtPlot': False,
     501                  'Yminmax': [Ymin, Ymax]}
     502    reader.Sample['ranId'] = valuesdict['ranId']
     503
     504    # Ending keys:
     505    # [u'Reflection Lists',
     506    #  u'Limits',
     507    #  'data',
     508    #  u'Index Peak List',
     509    #  u'Comments',
     510    #  u'Unit Cells List',
     511    #  u'Sample Parameters',
     512    #  u'Peak List',
     513    #  u'Background',
     514    #  u'Instrument Parameters']
     515    Tmin = np.min(reader.powderdata[0])
     516    Tmax = np.max(reader.powderdata[0])
     517
     518    default_background = [['chebyschev', False, 3, 1.0, 0.0, 0.0],
     519                          {'nDebye': 0, 'debyeTerms': [], 'nPeaks': 0, 'peaksList': []}]
     520
     521    output_dict = {u'Reflection Lists': {},
     522                   u'Limits': reader.pwdparms.get('Limits', [(Tmin, Tmax), [Tmin, Tmax]]),
     523                   u'data': [valuesdict, reader.powderdata, HistName],
     524                   u'Index Peak List': [[], []],
     525                   u'Comments': reader.comments,
     526                   u'Unit Cells List': [],
     527                   u'Sample Parameters': reader.Sample,
     528                   u'Peak List': {'peaks': [], 'sigDict': {}},
     529                   u'Background': reader.pwdparms.get('Background', default_background),
     530                   u'Instrument Parameters': [Iparm1, Iparm2],
     531                   }
     532
     533    names = [u'Comments',
     534             u'Limits',
     535             u'Background',
     536             u'Instrument Parameters',
     537             u'Sample Parameters',
     538             u'Peak List',
     539             u'Index Peak List',
     540             u'Unit Cells List',
     541             u'Reflection Lists']
     542
     543    # TODO controls?? GSASII.py:1664-7
     544
     545    return HistName, [HistName] + names, output_dict
     546
     547
     548def _deep_copy_into(from_, into):
     549    """Helper function for reloading .gpx file. See G2Project.reload()
     550    Author: Jackson O'Donnell jacksonhodonnell@gmail.com
     551    """
     552    if isinstance(from_, dict) and isinstance(into, dict):
     553        combined_keys = set(from_.keys()).union(into.keys())
     554        for key in combined_keys:
     555            if key in from_ and key in into:
     556                both_dicts = (isinstance(from_[key], dict)
     557                              and isinstance(into[key], dict))
     558                both_lists = (isinstance(from_[key], list)
     559                              and isinstance(into[key], list))
     560                if both_dicts or both_lists:
     561                    _deep_copy_into(from_[key], into[key])
     562                else:
     563                    into[key] = from_[key]
     564            elif key in from_:
     565                into[key] = from_[key]
     566            else:  # key in into
     567                del into[key]
     568    elif isinstance(from_, list) and isinstance(into, list):
     569        if len(from_) == len(into):
     570            for i in xrange(len(from_)):
     571                both_dicts = (isinstance(from_[i], dict)
     572                              and isinstance(into[i], dict))
     573                both_lists = (isinstance(from_[i], list)
     574                              and isinstance(into[i], list))
     575                if both_dicts or both_lists:
     576                    _deep_copy_into(from_[i], into[i])
     577                else:
     578                    into[i] = from_[i]
     579        else:
     580            into[:] = from_
     581
     582
     583class G2ObjectWrapper(object):
     584    """Base class for all GSAS-II object wrappers
     585    Author: Jackson O'Donnell jacksonhodonnell@gmail.com
     586    """
     587    def __init__(self, datadict):
     588        self.data = datadict
     589
     590    def __getitem__(self, key):
     591        return self.data[key]
     592
     593    def __setitem__(self, key, value):
     594        self.data[key] = value
     595
     596    def __contains__(self, key):
     597        return key in self.data
     598
     599    def get(self, k, d=None):
     600        return self.data.get(k, d)
     601
     602    def keys(self):
     603        return self.data.keys()
     604
     605    def values(self):
     606        return self.data.values()
     607
     608    def items(self):
     609        return self.data.items()
     610
     611
     612class G2Project(G2ObjectWrapper):
     613    def __init__(self, gpxfile=None, author=None, filename=None):
     614        """Loads a GSAS-II project from a specified filename.
     615
     616        :param str gpxfile: Existing .gpx file to be loaded. If nonexistent,
     617        creates an empty project.
     618        :param str author: Author's name. Optional.
     619        :param str filename: The filename the project should be saved to in
     620        the future. If both filename and gpxfile are present, the project is
     621        loaded from the gpxfile, then set to  save to filename in the future"""
     622        if gpxfile is None:
     623            self.filename = filename
     624            self.data, self.names = make_empty_project(author=author, filename=filename)
     625        elif isinstance(gpxfile, str):
     626            # TODO set author, filename
     627            self.filename = gpxfile
     628            self.data, self.names = LoadDictFromProjFile(gpxfile)
     629            self.index_ids()
     630        else:
     631            raise ValueError("Not sure what to do with gpxfile")
     632
     633    @classmethod
     634    def from_dict_and_names(cls, gpxdict, names, filename=None):
     635        out = cls()
     636        if filename:
     637            out.filename = filename
     638            gpxdict['Controls']['data']['LastSavedAs'] = filename
     639        else:
     640            try:
     641                out.filename = gpxdict['Controls']['data']['LastSavedAs']
     642            except KeyError:
     643                out.filename = None
     644        out.data = gpxdict
     645        out.names = names
     646
     647    def save(self, filename=None):
     648        """Saves the project, either to the current filename, or to a new file.
     649
     650        Updates self.filename if a new filename provided"""
     651        # TODO update LastSavedUsing ?
     652        if filename:
     653            self.data['Controls']['data']['LastSavedAs'] = filename
     654            self.filename = filename
     655        elif not self.filename:
     656            raise AttributeError("No file name to save to")
     657        SaveDictToProjFile(self.data, self.names, self.filename)
     658
     659    def add_powder_histogram(self, datafile, iparams):
     660        """Loads a powder data histogram into the project.
     661
     662        Automatically checks for an instrument parameter file, or one can be
     663        provided."""
     664        pwdrreaders = import_generic(datafile, PwdrDataReaders)
     665        histname, new_names, pwdrdata = load_pwd_from_reader(
     666                                          pwdrreaders[0], iparams,
     667                                          list(self.histogram_names()))
     668        if histname in self.data:
     669            print("Warning - redefining histogram", histname)
     670        else:
     671            if self.names[-1][0] == 'Phases':
     672                self.names.insert(-1, new_names)
     673            else:
     674                self.names.append(new_names)
     675        self.data[histname] = pwdrdata
     676        return self.histogram(histname)
     677
     678    def add_phase(self, phasefile, phasename=None, histograms=[]):
     679        """Loads a phase into the project from a .cif file
     680
     681        :returns: A :class:`gpx_manipulations.G2Phase` object representing the
     682        new phase."""
     683        # TODO handle multiple phases in a file
     684        phasereaders = import_generic(phasefile, PhaseReaders)
     685        phasereader = phasereaders[0]
     686
     687        phasename = phasename or phasereader.Phase['General']['Name']
     688        phaseNameList = list(self.phase_names())
     689        phasename = G2obj.MakeUniqueLabel(phasename, phaseNameList)
     690        phasereader.Phase['General']['Name'] = phasename
     691
     692        if 'Phases' not in self.data:
     693            self.data[u'Phases'] = { 'data': None }
     694        assert phasename not in self.data['Phases'], "phase names should be unique"
     695        self.data['Phases'][phasename] = phasereader.Phase
     696
     697        if phasereader.Constraints:
     698            Constraints = self.data['Constraints']
     699            for i in phasereader.Constraints:
     700                if isinstance(i, dict):
     701                    if '_Explain' not in Constraints:
     702                        Constraints['_Explain'] = {}
     703                    Constraints['_Explain'].update(i)
     704                else:
     705                    Constraints['Phase'].append(i)
     706
     707        data = self.data['Phases'][phasename]
     708        generalData = data['General']
     709        SGData = generalData['SGData']
     710        NShkl = len(G2spc.MustrainNames(SGData))
     711        NDij = len(G2spc.HStrainNames(SGData))
     712        Super = generalData.get('Super', 0)
     713        if Super:
     714            SuperVec = np.array(generalData['SuperVec'][0])
     715        else:
     716            SuperVec = []
     717        UseList = data['Histograms']
     718
     719        for histname in histograms:
     720            if histname.startswith('HKLF '):
     721                raise NotImplementedError("Does not support HKLF yet")
     722            elif histname.startswith('PWDR '):
     723                hist = self.histogram(histname)
     724                hist['Reflection Lists'][generalData['Name']] = {}
     725                UseList[histname] = SetDefaultDData('PWDR', histname, NShkl=NShkl, NDij=NDij)
     726                for key, val in [('Use', True), ('LeBail', False),
     727                                 ('newLeBail', True),
     728                                 ('Babinet', {'BabA': [0.0, False],
     729                                              'BabU': [0.0, False]})]:
     730                    if key not in UseList[histname]:
     731                        UseList[histname][key] = val
     732            else:
     733                raise NotImplementedError("Unexpected histogram" + histname)
     734
     735        for obj in self.names:
     736            if obj[0] == 'Phases':
     737                phasenames = obj
     738                break
     739        else:
     740            phasenames = [u'Phases']
     741            self.names.append(phasenames)
     742        phasenames.append(unicode(phasename))
     743
     744        # TODO should it be self.filename, not phasefile?
     745        SetupGeneral(data, os.path.dirname(phasefile))
     746        self.index_ids()
     747
     748        return self.phase(phasename)
     749
     750    def reload(self):
     751        """Reload self from self.filename"""
     752        data, names = LoadDictFromProjFile(self.filename)
     753        self.names = names
     754        # Need to deep copy the new data file data into the current tree,
     755        # so that any existing G2Phase, or G2PwdrData objects will still be
     756        # valid
     757        _deep_copy_into(from_=data, into=self.data)
     758
     759    def refine(self, newfile=None, printFile=None):
     760        # index_ids will automatically save the project
     761        self.index_ids()
     762        # TODO G2strMain does not properly use printFile
     763        G2strMain.Refine(self.filename)
     764        # Reload yourself
     765        self.reload()
     766
     767    def histogram_names(self):
     768        """Gives an iterator of the names of each histogram in the project."""
     769        for obj in self.names:
     770            if len(obj) > 1 and obj[0] != u'Phases':
     771                yield obj[0]
     772
     773    def histogram(self, histname):
     774        """Returns the histogram named histname, or None if it does not exist."""
     775        if histname in self.data:
     776            return G2PwdrData(self.data[histname])
     777        return None
     778
     779    def histograms(self):
     780        """Return a list of all histograms, as
     781        :class:`gpx_manipulations.G2PwdrData` objects"""
     782        return [self.histogram(name) for name in self.histogram_names()]
     783
     784    def phase_names(self):
     785        """Gives an iterator of the names of each phase in the project."""
     786        for obj in self.names:
     787            if obj[0] == 'Phases':
     788                for name in obj[1:]:
     789                    yield name
     790                break
     791
     792    def phase(self, phasename):
     793        return G2Phase(self.data['Phases'][phasename])
     794
     795    def phases(self):
     796        return [self.phase(name) for name in self.phase_names()]
     797
     798    def do_refinements(self, refinements, histogram='all', phase='all',
     799                       outputnames=None):
     800        """Conducts a series of refinements."""
     801        if outputnames:
     802            if len(refinements) != len(outputnames):
     803                raise ValueError("Should have same number of outuputs to"
     804                                 "refinements")
     805        else:
     806            outputnames = [None for r in refinements]
     807
     808        for output, refinement in zip(outputnames, refinements):
     809            self.set_refinement(refinement, histogram)
     810            # Handle 'once' args - refinements that are disabled after this
     811            # refinement
     812            if 'once' in refinement:
     813                temp = {'set': refinement['once']}
     814                self.set_refinement(temp, histogram, phase)
     815
     816            if output:
     817                self.save(output)
     818
     819            self.refine()  # newFile=output)
     820
     821            # Handle 'once' args - refinements that are disabled after this
     822            # refinement
     823            if 'once' in refinement:
     824                temp = {'clear': refinement['once']}
     825                self.set_refinement(temp, histogram, phase)
     826
     827    def set_refinement(self, refinement, histogram='all', phase='all'):
     828        """Apply specified refinements to a given histogram(s) or phase(s).
     829
     830        :param dict refinement: The refinements to be conducted
     831        :param histogram: Either a name of a histogram (str), a list of
     832        histogram names, or 'all' (default)
     833        :param phase: Either a name of a phase (str), a list of phase names, or
     834        'all' (default)"""
     835        if histogram == 'all':
     836            hists = [self.histogram(name)
     837                     for name in self.histogram_names()]
     838        elif isinstance(histogram, str):
     839            hists = [self.histogram(histogram)]
     840        else:
     841            hists = [self.histogram(name) for name in histogram]
     842
     843        if phase == 'all':
     844            phases = [self.phase(name) for name in self.phase_names()]
     845        elif isinstance(phase, str):
     846            phases = [self.phase(phase)]
     847        else:
     848            phases = [self.phase(name) for name in phase]
     849
     850        for key, val in refinement.get('set', {}).items():
     851            # Apply refinement options
     852            if G2PwdrData.is_valid_refinement_key(key):
     853                for hist in hists:
     854                    hist.set_refinement(key, val)
     855            elif G2Phase.is_valid_refinement_key(key):
     856                for phase in phases:
     857                    phase.set_refinement(key, val)
     858            else:
     859                print("Unknown refinement key:", key)
     860
     861        for key, val in refinement.get('clear', {}).items():
     862            # Clear refinement options
     863            if G2PwdrData.is_valid_refinement_key(key):
     864                for hist in hists:
     865                    hist.clear_refinement(key, val)
     866            elif G2Phase.is_valid_refinement_key(key):
     867                for phase in phases:
     868                    phase.clear_refinement(key, val)
     869            else:
     870                print("Unknown refinement key:", key)
     871
     872    def index_ids(self):
     873        self.save()
     874        return G2strIO.GetUsedHistogramsAndPhases(self.filename)
     875
     876    def add_constraint_raw(self, cons_scope, constr):
     877        """Adds a constraint of type consType to the project.
     878        cons_scope should be one of "Hist", "Phase", "HAP", or "Global".
     879
     880        WARNING it does not check the constraint is well-constructed"""
     881        constrs = self['Constraints']['data']
     882        if 'Global' not in constrs:
     883            constrs['Global'] = []
     884        constrs[cons_scope].append(constr)
     885
     886    def hold_many(self, vars, type):
     887        """Apply holds for all the variables in vars, for constraint of a given type.
     888
     889        type is passed directly to add_constraint_raw as consType
     890
     891        :param list vars: A list of variables to hold. Either G2obj.G2VarObj objects,
     892        string variable specifiers, or arguments for :meth:`make_var_obj`
     893        :param str type: A string constraint type specifier. See
     894        :meth:`add_constraint_raw`"""
     895        for var in vars:
     896            if isinstance(var, str):
     897                var = self.make_var_obj(var)
     898            elif not isinstance(var, G2obj.G2VarObj):
     899                var = self.make_var_obj(*var)
     900            self.add_constraint_raw(type, [[1.0, var], None, None, 'h'])
     901
     902    def make_var_obj(self, phase=None, hist=None, varname=None, atomId=None,
     903                     reloadIdx=True):
     904        """Wrapper to create a G2VarObj. Takes either a string representaiton ("p:h:name:a")
     905        or individual names of phase, histogram, varname, and atomId.
     906
     907        Automatically converts string phase, hist, or atom names into the ID required
     908        by G2VarObj."""
     909        if reloadIdx:
     910            self.index_ids()
     911
     912        # If string representation, short circuit
     913        if hist is None and varname is None and atomId is None:
     914            if isinstance(phase, str) and ':' in phase:
     915                return G2obj.G2VarObj(phase)
     916
     917        # Get phase index
     918        phaseObj = None
     919        if isinstance(phase, G2Phase):
     920            phaseObj = phase
     921            phase = G2obj.PhaseRanIdLookup[phase.ranId]
     922        elif phase in self['Phases']:
     923            phaseObj = self.phase(phase)
     924            phaseRanId = phaseObj.ranId
     925            phase = G2obj.PhaseRanIdLookup[phaseRanId]
     926
     927        # Get histogram index
     928        if isinstance(hist, G2PwdrData):
     929            hist = G2obj.HistRanIdLookup[hist.ranId]
     930        elif hist in self.data:
     931            histRanId = self.histogram(hist).ranId
     932            hist = G2obj.HistRanIdLookup[histRanId]
     933
     934        # Get atom index (if any)
     935        if isinstance(atomId, G2AtomRecord):
     936            atomId = G2obj.AtomRanIdLookup[phase][atomId.ranId]
     937        elif phaseObj:
     938            atomObj = phaseObj.atom(atomId)
     939            if atomObj:
     940                atomRanId = atomObj.ranId
     941                atomId = G2obj.AtomRanIdLookup[phase][atomRanId]
     942
     943        return G2obj.G2VarObj(phase, hist, varname, atomId)
     944
     945
     946class G2AtomRecord(G2ObjectWrapper):
     947    """Wrapper for an atom record
     948    Author: Jackson O'Donnell jacksonhodonnell@gmail.com
     949    """
     950    def __init__(self, data, indices):
     951        self.data = data
     952        self.cx, self.ct, self.cs, self.cia = indices
     953
     954    @property
     955    def label(self):
     956        return self.data[self.ct-1]
     957
     958    @property
     959    def type(self):
     960        return self.data[self.ct]
     961
     962    @property
     963    def refinement_flags(self):
     964        return self.data[self.ct+1]
     965
     966    @refinement_flags.setter
     967    def set_refinement(self, other):
     968        for c in other:
     969            if c not in ' FXU':
     970                raise ValueError("Invalid atom refinement: ", other)
     971        self.data[self.ct+1] = unicode(other)
     972
     973    @property
     974    def coordinates(self):
     975        return tuple(self.data[self.cx:self.cx+3])
     976
     977    @property
     978    def ranId(self):
     979        return self.data[self.cia+8]
     980
     981    @property
     982    def adp_flag(self):
     983        # Either 'I' or 'A'
     984        return self.data[self.cia]
     985
     986    @property
     987    def uiso(self):
     988        if self.adp_flag == 'I':
     989            return self.data[self.cia+1]
     990        else:
     991            return self.data[self.cia+2:self.cia+8]
     992
     993
     994class G2PwdrData(G2ObjectWrapper):
     995    """Wraps a Poweder Data Histogram.
     996    Author: Jackson O'Donnell jacksonhodonnell@gmail.com
     997    """
     998    @staticmethod
     999    def is_valid_refinement_key(key):
     1000        valid_keys = ['Limits', 'Sample Parameters', 'Background',
     1001                      'Instrument Parameters']
     1002        return key in valid_keys
     1003
     1004    @property
     1005    def name(self):
     1006        return self['data'][-1]
     1007
     1008    @property
     1009    def ranId(self):
     1010        return self['data'][0]['ranId']
     1011
     1012    def fit_fixed_points(self):
     1013        """Attempts to apply a background fit to the fixed points currently specified."""
     1014
     1015        def SetInstParms(Inst):
     1016            dataType = Inst['Type'][0]
     1017            insVary = []
     1018            insNames = []
     1019            insVals = []
     1020            for parm in Inst:
     1021                insNames.append(parm)
     1022                insVals.append(Inst[parm][1])
     1023                if parm in ['U','V','W','X','Y','SH/L','I(L2)/I(L1)','alpha',
     1024                    'beta-0','beta-1','beta-q','sig-0','sig-1','sig-2','sig-q',] and Inst[parm][2]:
     1025                        Inst[parm][2] = False
     1026            instDict = dict(zip(insNames, insVals))
     1027            instDict['X'] = max(instDict['X'], 0.01)
     1028            instDict['Y'] = max(instDict['Y'], 0.01)
     1029            if 'SH/L' in instDict:
     1030                instDict['SH/L'] = max(instDict['SH/L'], 0.002)
     1031            return dataType, instDict, insVary
     1032
     1033        bgrnd = self['Background']
     1034
     1035        # Need our fixed points in order
     1036        bgrnd[1]['FixedPoints'].sort(key=lambda pair: pair[0])
     1037        X = [x for x, y in bgrnd[1]['FixedPoints']]
     1038        Y = [y for x, y in bgrnd[1]['FixedPoints']]
     1039
     1040        # Some simple lookups
     1041        limits = self['Limits']
     1042        inst, inst2 = self['Instrument Parameters']
     1043        pwddata = self['data'][1]
     1044
     1045        # Construct the data for background fitting
     1046        xBeg = np.searchsorted(pwddata[0], limits[0])
     1047        xFin = np.searchsorted(pwddata[0], limits[1])
     1048        xdata = pwddata[0][xBeg:xFin]
     1049        ydata = si.interp1d(X,Y)(ma.getdata(xdata))
     1050
     1051        W = [1]*len(xdata)
     1052        Z = [0]*len(xdata)
     1053
     1054        dataType, insDict, insVary = SetInstParms(inst)
     1055        bakType, bakDict, bakVary = G2pwd.SetBackgroundParms(bgrnd)
     1056
     1057        # Do the fit
     1058        # TODO get controls from root
     1059        data = np.array([xdata, ydata, W, Z, Z, Z])
     1060        G2pwd.DoPeakFit('LSQ', [], bgrnd, limits, inst, inst2, data,
     1061                        preVaryList=bakVary, controls={})
     1062
     1063        # Post-fit
     1064        parmDict = {}
     1065        bakType, bakDict, bakVary = G2pwd.SetBackgroundParms(bgrnd)
     1066        parmDict.update(bakDict)
     1067        parmDict.update(insDict)
     1068        pwddata[3][xBeg:xFin] *= 0
     1069        pwddata[5][xBeg:xFin] *= 0
     1070        pwddata[4][xBeg:xFin] = G2pwd.getBackground('', parmDict, bakType, dataType, xdata)[0]
     1071
     1072        # TODO adjust pwddata? GSASIIpwdGUI.py:1041
     1073
     1074    def y_calc(self):
     1075        return self['data'][1][3]
     1076
     1077    def set_refinement(self, key, value):
     1078        """Sets the refinement parameter 'key' to the specification 'value'"""
     1079        if key == 'Limits':
     1080            old_limits = self.data['Limits'][1]
     1081            new_limits = value
     1082            if isinstance(new_limits, dict):
     1083                if 'low' in new_limits:
     1084                    old_limits[0] = new_limits['low']
     1085                if 'high' in new_limits:
     1086                    old_limits[1] = new_limits['high']
     1087            else:
     1088                old_limits[0], old_limits[1] = new_limits
     1089        elif key == 'Sample Parameters':
     1090            sample = self.data['Sample Parameters']
     1091            for sparam in value:
     1092                sample[sparam][1] = True
     1093        elif key == 'Background':
     1094            bkg, peaks = self.data['Background']
     1095
     1096            # If True or False, just set the refine parameter
     1097            if value in (True, False):
     1098                bkg[1] = value
     1099                return
     1100
     1101            if 'type' in value:
     1102                bkg[0] = value['type']
     1103            if 'refine' in value:
     1104                bkg[1] = value['refine']
     1105            if 'no. coeffs' in value:
     1106                cur_coeffs = bkg[2]
     1107                n_coeffs = value['no. coeffs']
     1108                if n_coeffs > cur_coeffs:
     1109                    for x in range(n_coeffs - cur_coeffs):
     1110                        bkg.append(0.0)
     1111                else:
     1112                    for _ in range(cur_coeffs - n_coeffs):
     1113                        bkg.pop()
     1114                bkg[2] = n_coeffs
     1115            if 'coeffs' in value:
     1116                bkg[3:] = value['coeffs']
     1117            if 'FixedPoints' in value:
     1118                peaks['FixedPoints'] = value
     1119            if value.get('fit fixed points', False):
     1120                self.fit_fixed_points()
     1121
     1122        elif key == 'Instrument Parameters':
     1123            instrument, secondary = self.data['Instrument Parameters']
     1124            for iparam in value:
     1125                try:
     1126                    instrument[iparam][2] = True
     1127                except IndexError:
     1128                    raise ValueError("Invalid key:", iparam)
     1129        else:
     1130            raise ValueError("Unknown key:", key)
     1131
     1132    def clear_refinement(self, key, value):
     1133        """Clears the refinement parameter 'key' and its associated value."""
     1134        if key == 'Limits':
     1135            old_limits, cur_limits = self.data['Limits']
     1136            cur_limits[0], cur_limits[1] = old_limits
     1137        elif key == 'Sample Parameters':
     1138            sample = self.data['Sample Parameters']
     1139            for sparam in value:
     1140                sample[sparam][1] = False
     1141        elif key == 'Background':
     1142            bkg, peaks = self.data['Background']
     1143
     1144            # If True or False, just set the refine parameter
     1145            if value in (True, False):
     1146                bkg[1] = False
     1147                return
     1148
     1149            bkg[1] = False
     1150            if 'FixedPoints' in value:
     1151                if 'FixedPoints' in peaks:
     1152                    del peaks['FixedPoints']
     1153        elif key == 'Instrument Parameters':
     1154            instrument, secondary = self.data['Instrument Parameters']
     1155            for iparam in value:
     1156                instrument[iparam][2] = False
     1157        else:
     1158            raise ValueError("Unknown key:", key)
     1159
     1160
     1161class G2Phase(G2ObjectWrapper):
     1162    """A wrapper object around a given phase.
     1163    Author: Jackson O'Donnell jacksonhodonnell@gmail.com
     1164    """
     1165    @staticmethod
     1166    def is_valid_refinement_key(key):
     1167        valid_keys = ["Cell", "Atoms", "LeBail"]
     1168        return key in valid_keys
     1169
     1170    def atom(self, atomlabel):
     1171        """Returns the atom specified by atomlabel, or None if it does not
     1172        exist."""
     1173        # Consult GSASIIobj.py for the meaning of this
     1174        cx, ct, cs, cia = self.data['General']['AtomPtrs']
     1175        ptrs = [cx, ct, cs, cia]
     1176        atoms = self.data['Atoms']
     1177        for atom in atoms:
     1178            if atom[ct-1] == atomlabel:
     1179                return G2AtomRecord(atom, ptrs)
     1180
     1181    def atoms(self):
     1182        """Returns a list of atoms present in the phase."""
     1183        ptrs = self.data['General']['AtomPtrs']
     1184        output = []
     1185        atoms = self.data['Atoms']
     1186        for atom in atoms:
     1187            output.append(G2AtomRecord(atom, ptrs))
     1188        return output
     1189
     1190    @property
     1191    def ranId(self):
     1192        return self.data['ranId']
     1193
     1194    def cell_dict(self):
     1195        cell = self['General']['Cell']
     1196        return {'a': cell[1], 'b': cell[2], 'c': cell[3],
     1197                'alpha': cell[4], 'beta': cell[5], 'gamma': cell[6],
     1198                'vol': cell[6]}
     1199
     1200    def set_refinement(self, key, value):
     1201        if key == "Cell":
     1202            self.data['General']['Cell'][0] = True
     1203        elif key == "Atoms":
     1204            cx, ct, cs, cia = self.data['General']['AtomPtrs']
     1205
     1206            for atomlabel, atomrefinement in value.items():
     1207                atom = self.atom(atomlabel)
     1208                atom.refinement_flags = atomrefinement
     1209        elif key == "LeBail":
     1210            hists = self.data['Histograms']
     1211            for hname, hoptions in hists.items():
     1212                if 'LeBail' not in hoptions:
     1213                    hoptions['newLeBail'] = True
     1214                hoptions['LeBail'] = bool(value)
     1215        else:
     1216            raise ValueError("Unknown key:", key)
     1217
     1218    def clear_refinement(self, key, value):
     1219        if key == "Cell":
     1220            self.data['General']['Cell'][0] = False
     1221        elif key == "Atoms":
     1222            cx, ct, cs, cia = self.data['General']['AtomPtrs']
     1223
     1224            for atomlabel in value:
     1225                atom = self.atom(atomlabel)
     1226                # Set refinement to none
     1227                atom.refinement_flags = ' '
     1228        elif key == "LeBail":
     1229            hists = self.data['Histograms']
     1230            for hname, hoptions in hists.items():
     1231                if 'LeBail' not in hoptions:
     1232                    hoptions['newLeBail'] = True
     1233                hoptions['LeBail'] = False
     1234        else:
     1235            raise ValueError("Unknown key:", key)
    1421236
    1431237def main():
    144     'Needs a doc string'
     1238    '''TODO: the command-line options need some thought
     1239    '''
    1451240    arg = sys.argv
    146     print arg
     1241    print(arg)
    1471242    if len(arg) > 1:
    1481243        GPXfile = arg[1]
    1491244        if not ospath.exists(GPXfile):
    150             print 'ERROR - ',GPXfile," doesn't exist!"
     1245            print(u'ERROR - '+GPXfile+u" doesn't exist!")
    1511246            exit()
    1521247        Project,nameList = LoadDictFromProjFile(GPXfile)
    1531248        SaveDictToProjFile(Project,nameList,'testout.gpx')
    1541249    else:
    155         print 'ERROR - missing filename'
     1250        print('ERROR - missing filename')
    1561251        exit()
    157     print "Done"
     1252    print("Done")
    1581253         
    1591254if __name__ == '__main__':
     1255    PwdrDataReaders = G2fil.LoadImportRoutines("pwd", "Powder_Data")
     1256    PhaseReaders = G2fil.LoadImportRoutines("phase", "Phase")
    1601257    main()
     1258
     1259    # from gpx_manipulatons.py
     1260    # try:
     1261    #     filename, authorname = sys.argv[1:3]
     1262    #     proj, names = make_empty_project(authorname, filename)
     1263    #     SaveDictToProjFile(proj, names, os.path.abspath(filename))
     1264    # except ValueError:
     1265    #     print("Usage: {} <filename> <author>".format(sys.argv[0]))
     1266    #     sys.exit(-1)
     1267
     1268
     1269    # from refinements.py
     1270#      USAGE = """USAGE: {} datafile instparams phasefile projectname refinements
     1271
     1272# datafile:    Input powder data
     1273# intparams:   Corresponding instrument parameter file
     1274# phasefile:   Phase to refine against data
     1275# projectname: Project file to be created, should end in .gpx
     1276# refinements: JSON file of refinements to be executed
     1277# """
     1278#     try:
     1279#         datafile, instprm, phasefile, projectname, refinements = sys.argv[1:]
     1280#     except ValueError:
     1281#         print(USAGE.format(sys.argv[0]))
     1282#         sys.exit(-1)
     1283
     1284#     try:
     1285#         with open(refinements) as f:
     1286#             refinements = json.load(f)
     1287#     except IOError:
     1288#         print("No such refinements file: {}".format(refinements))
     1289
     1290#     print("Creating project file \"{}\"...".format(projectname))
     1291#     proj = G2Project(filename=projectname)
     1292#     # Add the histogram
     1293#     hist = proj.add_powder_histogram(datafile, instprm)
     1294#     # Add the phase, and associate it with the histogram
     1295#     proj.add_phase(phasefile, histograms=[hist.name])
     1296
     1297#     proj.do_refinements(refinements['refinements'])
     1298#     proj.save()
     1299
     1300
     1301    # from gpx_dumper
     1302    # import IPython.lib.pretty as pretty
     1303    # proj, nameList = LoadDictFromProjFile(sys.argv[1])
     1304    # print("names:", nameList)
     1305    # for key, val in proj.items():
     1306    #     print(key, ":", sep='')
     1307    #     pretty.pprint(val)
     1308   
Note: See TracChangeset for help on using the changeset viewer.