Changeset 2681 for trunk/GSASIIpwd.py


Ignore:
Timestamp:
Jan 31, 2017 3:38:48 PM (5 years ago)
Author:
vondreele
Message:

fix missing 'Oblique' problem in G2imgGUI
fix inner 2-theta limit issue in OnTransferAngles?
PDF Peaks can now clear all peaks
PDF peak picking now gives pos,mag & sig=0.085 as default
calculated PDF peak fit curve now plotted on G(R) from PDF Peaks
fitting of PDF peaks implemented - works fine; needs text output

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/GSASIIpwd.py

    r2660 r2681  
    1818import os
    1919import subprocess as subp
     20import copy
    2021
    2122import numpy as np
     
    286287    '''
    287288    auxPlot = []
    288     import copy
    289289    import scipy.fftpack as ft
    290290    Ibeg = np.searchsorted(xydata['Sample'][1][0],limits[0])
     
    373373        xydata['GofR'][1][1] = np.where(xydata['GofR'][1][0]<0.5,0.,xydata['GofR'][1][1])
    374374    return auxPlot
     375   
     376def PDFPeakFit(peaks,data):
     377    rs2pi = 1./np.sqrt(2*np.pi)
     378   
     379    def MakeParms(peaks):
     380        varyList = []
     381        parmDict = {'slope':peaks['Background'][1][1]}
     382        if peaks['Background'][2]:
     383            varyList.append('slope')
     384        for i,peak in enumerate(peaks['Peaks']):
     385            parmDict[str(i)+':pos'] = peak[0]
     386            parmDict[str(i)+':mag'] = peak[1]
     387            parmDict[str(i)+':sig'] = peak[2]
     388            if 'P' in peak[3]:
     389                varyList.append(str(i)+':pos')
     390            if 'M' in peak[3]:
     391                varyList.append(str(i)+':mag')
     392            if 'S' in peak[3]:
     393                varyList.append(str(i)+':sig')
     394        return parmDict,varyList
     395       
     396    def SetParms(peaks,parmDict,varyList):
     397        if 'slope' in varyList:
     398            peaks['Background'][1][1] = parmDict['slope']
     399        for i,peak in enumerate(peaks['Peaks']):
     400            if str(i)+':pos' in varyList:
     401                peak[0] = parmDict[str(i)+':pos']
     402            if str(i)+':mag' in varyList:
     403                peak[1] = parmDict[str(i)+':mag']
     404            if str(i)+':sig' in varyList:
     405                peak[2] = parmDict[str(i)+':sig']
     406       
     407   
     408    def CalcPDFpeaks(parmdict,Xdata):
     409        Z = parmDict['slope']*Xdata
     410        ipeak = 0
     411        while True:
     412            try:
     413                pos = parmdict[str(ipeak)+':pos']
     414                mag = parmdict[str(ipeak)+':mag']
     415                wid = parmdict[str(ipeak)+':sig']
     416                wid2 = 2.*wid**2
     417                Z += mag*rs2pi*np.exp(-(Xdata-pos)**2/wid2)/wid
     418                ipeak += 1
     419            except KeyError:        #no more peaks to process
     420                return Z
     421               
     422    def errPDFProfile(values,xdata,ydata,parmdict,varylist):       
     423        parmdict.update(zip(varylist,values))
     424        M = CalcPDFpeaks(parmdict,xdata)-ydata
     425        return M
     426                               
     427           
     428    print 'Do PDF peak fitting'
     429    newpeaks = copy.copy(peaks)
     430    iBeg = np.searchsorted(data[1][0],newpeaks['Limits'][0])
     431    iFin = np.searchsorted(data[1][0],newpeaks['Limits'][1])
     432    X = data[1][0][iBeg:iFin]
     433    Y = data[1][1][iBeg:iFin]
     434    parmDict,varyList = MakeParms(peaks)
     435    if not len(varyList):
     436        print ' Nothing varied'
     437        return newpeaks
     438   
     439    begin = time.time()
     440    values =  np.array(Dict2Values(parmDict, varyList))
     441    result = so.leastsq(errPDFProfile,values,full_output=True,ftol=0.0001,
     442           args=(X,Y,parmDict,varyList))
     443    ncyc = int(result[2]['nfev']/2)
     444    runtime = time.time()-begin   
     445    print 'PDF fitpeak time = %8.3fs, %8.3fs/cycle'%(runtime,runtime/ncyc)
     446#    chisq = np.sum(result[2]['fvec']**2)
     447    Values2Dict(parmDict, varyList, result[0])
     448    SetParms(peaks,parmDict,varyList)
     449   
     450    Z = CalcPDFpeaks(parmDict,X)
     451    newpeaks['calc'] = [X,Z]
     452    return newpeaks   
    375453   
    376454def MakeRDF(RDFcontrols,background,inst,pwddata):
Note: See TracChangeset for help on using the changeset viewer.