Changeset 686


Ignore:
Timestamp:
Jul 11, 2012 4:01:20 PM (9 years ago)
Author:
vondreele
Message:

add restraints to GSASIIstruct.py - input & printing only
no calcs yet

Location:
trunk
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • trunk/GSASII.py

    r683 r686  
    20392039        Histograms,Phases = self.GetUsedHistogramsAndPhasesfromTree()
    20402040        print Histograms.keys()
    2041         Natoms,phaseVary,phaseDict,pawleyLookup,FFtable,BLtable = G2str.GetPhaseData(Phases,Print=False)       
     2041        Natoms,phaseVary,phaseDict,pawleyLookup,FFtable,BLtable = G2str.GetPhaseData(Phases,RestDict=None,Print=False)       
    20422042        hapVary,hapDict,controlDict = G2str.GetHistogramPhaseData(Phases,Histograms,Print=False)
    20432043        histVary,histDict,controlDict = G2str.GetHistogramData(Histograms,Print=False)
  • trunk/GSASIIstruct.py

    r679 r686  
    128128    if not Histograms:
    129129        return 'Error: no diffraction data',''
    130     Natoms,phaseVary,phaseDict,pawleyLookup,FFtables,BLtables = GetPhaseData(Phases,Print=False)
     130    Natoms,phaseVary,phaseDict,pawleyLookup,FFtables,BLtables = GetPhaseData(Phases,RestDict=None,Print=False)
    131131    hapVary,hapDict,controlDict = GetHistogramPhaseData(Phases,Histograms,Print=False)
    132132    histVary,histDict,controlDict = GetHistogramData(Histograms,Print=False)
     
    134134    constrDict,fixedList = GetConstraints(GPXfile)
    135135    return G2mv.CheckConstraints(varyList,constrDict,fixedList)
     136   
     137def GetRestraints(GPXfile):
     138    '''Read the restraints from the GPX file
     139    '''
     140    fl = open(GPXfile,'rb')
     141    while True:
     142        try:
     143            data = cPickle.load(fl)
     144        except EOFError:
     145            break
     146        datum = data[0]
     147        if datum[0] == 'Restraints':
     148            restDict = datum[1]
     149    fl.close()
     150    return restDict
    136151   
    137152def GetPhaseNames(GPXfile):
     
    560575################################################################################       
    561576                   
    562 def GetPhaseData(PhaseData,Print=True):
     577def GetPhaseData(PhaseData,RestDict=None,Print=True):
    563578           
    564579    def PrintFFtable(FFtable):
     
    633648            ptstr += '%12.4f'%(SHcoeff[item])
    634649        print ptlbls
    635         print ptstr   
     650        print ptstr
     651       
     652    def PrintRestraints(phaseRest):
     653        if phaseRest:
     654            print '\n Restraints:'
     655            names = ['Bonds','Angles','Planes','Volumes']
     656            for i,name in enumerate(['Bond','Angle','Plane','Chiral']):
     657                itemRest = phaseRest[name]
     658                if itemRest[names[i]]:
     659                    print '\n  %30s %10.3f'%(name+' restraint weight factor',itemRest['wtFactor'])
     660                    print '  atoms(symOp), calc, obs, sig: '
     661                    for item in phaseRest[name][names[i]]:
     662                        text = '   '
     663                        for a,at in enumerate(item[0]):
     664                            text += at+'+('+item[1][a]+') '
     665                            if (a+1)%5 == 0:
     666                                text += '\n   '
     667                        print text,' %.3f %.3f %.3f'%(item[3],item[4],item[5])
    636668       
    637669    if Print: print ' Phases:'
     
    741773                    '%.3f'%(cell[6]),' volume =','%.3f'%(cell[7]),' Refine?',cell[0]
    742774                PrintTexture(textureData)
     775                if name in RestDict:
     776                    PrintRestraints(RestDict[name])
    743777                   
    744778        elif PawleyRef:
     
    29152949    parmdict.update(zip(varylist,values))
    29162950    G2mv.Dict2Map(parmdict,varylist)
    2917     Histograms,Phases = HistoPhases
     2951    Histograms,Phases,restDict = HistoPhases
    29182952    nvar = len(varylist)
    29192953    dMdv = np.empty(0)
     
    29803014    parmdict.update(zip(varylist,values))
    29813015    G2mv.Dict2Map(parmdict,varylist)
    2982     Histograms,Phases = HistoPhases
     3016    Histograms,Phases,restDict = HistoPhases
    29833017    nvar = len(varylist)
    29843018    Hess = np.empty(0)
     
    30663100    Values2Dict(parmdict, varylist, values)
    30673101    G2mv.Dict2Map(parmdict,varylist)
    3068     Histograms,Phases = HistoPhases
     3102    Histograms,Phases,restDict = HistoPhases
    30693103    M = np.empty(0)
    30703104    SumwYo = 0
     
    31843218    calcControls.update(Controls)           
    31853219    constrDict,fixedList = GetConstraints(GPXfile)
     3220    restDict = GetRestraints(GPXfile)
    31863221    Histograms,Phases = GetUsedHistogramsAndPhases(GPXfile)
    31873222    if not Phases:
     
    31933228        print ' *** Refine aborted ***'
    31943229        raise Exception       
    3195     Natoms,phaseVary,phaseDict,pawleyLookup,FFtables,BLtables = GetPhaseData(Phases)
     3230    Natoms,phaseVary,phaseDict,pawleyLookup,FFtables,BLtables = GetPhaseData(Phases,restDict)
    31963231    calcControls['Natoms'] = Natoms
    31973232    calcControls['FFtables'] = FFtables
     
    32373272            result = so.leastsq(errRefine,values,Dfun=dervRefine,full_output=True,
    32383273                ftol=Ftol,col_deriv=True,factor=Factor,
    3239                 args=([Histograms,Phases],parmDict,varyList,calcControls,pawleyLookup,dlg))
     3274                args=([Histograms,Phases,restDict],parmDict,varyList,calcControls,pawleyLookup,dlg))
    32403275            ncyc = int(result[2]['nfev']/2)
    32413276        elif 'Hessian' in Controls['deriv type']:
    32423277            result = G2mth.HessianLSQ(errRefine,values,Hess=HessRefine,ftol=Ftol,maxcyc=maxCyc,
    3243                 args=([Histograms,Phases],parmDict,varyList,calcControls,pawleyLookup,dlg))
     3278                args=([Histograms,Phases,restDict],parmDict,varyList,calcControls,pawleyLookup,dlg))
    32443279            ncyc = result[2]['num cyc']+1
    32453280            Rvals['lamMax'] = result[2]['lamMax']
    32463281        else:           #'numeric'
    32473282            result = so.leastsq(errRefine,values,full_output=True,ftol=Ftol,epsfcn=1.e-8,factor=Factor,
    3248                 args=([Histograms,Phases],parmDict,varyList,calcControls,pawleyLookup,dlg))
     3283                args=([Histograms,Phases,restDict],parmDict,varyList,calcControls,pawleyLookup,dlg))
    32493284            ncyc = int(result[2]['nfev']/len(varyList))
    32503285#        table = dict(zip(varyList,zip(values,result[0],(result[0]-values))))
     
    33363371    ShowControls(Controls)           
    33373372    constrDict,fixedList = GetConstraints(GPXfile)
     3373    restDict = GetRestraints(GPXfile)
    33383374    Histograms,Phases = GetUsedHistogramsAndPhases(GPXfile)
    33393375    if not Phases:
     
    33453381        print ' *** Refine aborted ***'
    33463382        raise Exception
    3347     Natoms,phaseVary,phaseDict,pawleyLookup,FFtables,BLtables = GetPhaseData(Phases,False)
     3383    Natoms,phaseVary,phaseDict,pawleyLookup,FFtables,BLtables = GetPhaseData(Phases,restDict,False)
    33483384    if 'Seq Data' in Controls:
    33493385        histNames = Controls['Seq Data']
     
    34233459                result = so.leastsq(errRefine,values,Dfun=dervRefine,full_output=True,
    34243460                    ftol=Ftol,col_deriv=True,factor=Factor,
    3425                     args=([Histo,Phases],parmDict,varyList,calcControls,pawleyLookup,dlg))
     3461                    args=([Histo,Phases,restDict],parmDict,varyList,calcControls,pawleyLookup,dlg))
    34263462                ncyc = int(result[2]['nfev']/2)
    34273463            elif 'Hessian' in Controls['deriv type']:
    34283464                result = G2mth.HessianLSQ(errRefine,values,Hess=HessRefine,ftol=Ftol,maxcyc=maxCyc,
    3429                     args=([Histo,Phases],parmDict,varyList,calcControls,pawleyLookup,dlg))
     3465                    args=([Histo,Phases,restDict],parmDict,varyList,calcControls,pawleyLookup,dlg))
    34303466                ncyc = result[2]['num cyc']+1                           
    34313467            else:           #'numeric'
    34323468                result = so.leastsq(errRefine,values,full_output=True,ftol=Ftol,epsfcn=1.e-8,factor=Factor,
    3433                     args=([Histo,Phases],parmDict,varyList,calcControls,pawleyLookup,dlg))
     3469                    args=([Histo,Phases,restDict],parmDict,varyList,calcControls,pawleyLookup,dlg))
    34343470                ncyc = int(result[2]['nfev']/len(varyList))
    34353471
Note: See TracChangeset for help on using the changeset viewer.