Changeset 2129 for trunk


Ignore:
Timestamp:
Jan 21, 2016 10:25:48 AM (9 years ago)
Author:
vondreele
Message:

revisions to FWHM calculations - some math errors fixed
implement merge fully - Laue groups still need work
triclinic, monoclinic & orthorhombic OK, others not

Location:
trunk
Files:
5 edited

Legend:

Unmodified
Added
Removed
  • TabularUnified trunk/GSASII.py

    r2127 r2129  
    705705                                    ref[3+Super] = 0
    706706                                else:
    707                                     ref[3+Super] = 1    #twin id?
     707                                    ref[3+Super] = 1    #twin id
    708708                            else:
    709709                                ref[3] = 0
     
    31663166                            dsp = G2lat.Pos2dsp(Inst,peak[0])
    31673167                            if 'T' in Type:  #TOF - more cols
    3168                                 FWHM = 2.*G2pwd.getgamFW(peak[10],peak[8])      #to get delta-TOF from Gam(peak)
     3168                                sig = 2.35482*np.sqrt(peak[8])
     3169                                gam = peak[10]
     3170                                FWHM = 2.*G2pwd.getgamFW(gam,sig)      #to get delta-TOF from Gam(peak)
    31693171                                file.write("%10.2f %10.5f %12.2f %10.3f %10.3f %10.3f %10.3f %10.3f\n" % \
    31703172                                    (peak[0],dsp,peak[2],np.sqrt(max(0.0001,peak[4])),peak[6],peak[8],peak[10],FWHM))
    31713173                            else:               #CW
    3172                                 FWHM = 2.*G2pwd.getgamFW(peak[6],peak[4])      #to get delta-2-theta in deg. from Gam(peak)
     3174                                sig = 2.35482*np.sqrt(peak[4])/100.
     3175                                gam = peak[6]/100.
     3176                                FWHM = 2.*G2pwd.getgamFW(gam,sig)      #to get delta-2-theta in deg. from Gam(peak)
    31733177                                file.write("%10.3f %10.5f %12.2f %10.5f %10.5f %10.5f \n" % \
    31743178                                    (peak[0],dsp,peak[2],np.sqrt(max(0.0001,peak[4]))/100.,peak[6]/100.,FWHM/100.)) #convert to deg
     
    32053209                                        file.write('%s \n'%('   h   k   l   m    d-space   2-theta       wid        F**2'))
    32063210                                    for peak in peaks['RefList']:
    3207                                         FWHM = 2.*G2pwd.getgamFW(peak[7],peak[6])
    32083211                                        if 'T' in peaks.get('Type','PXC'):
     3212                                            sig = 2.35482*np.sqrt(peak[6])
     3213                                            gam = peak[7]
     3214                                            FWHM = 2.*G2pwd.getgamFW(gam,sig)
    32093215                                            file.write(" %3d %3d %3d %3d %10.5f %10.2f %10.5f %10.3f \n" % \
    32103216                                                (int(peak[0]),int(peak[1]),int(peak[2]),int(peak[3]),peak[4],peak[5],FWHM,peak[8]))
    32113217                                        else:
     3218                                            sig = 2.35482*np.sqrt(peak[6])/100.
     3219                                            gam = peak[7]/100.
     3220                                            FWHM = 2.*G2pwd.getgamFW(gam,sig)
    32123221                                            file.write(" %3d %3d %3d %3d %10.5f %10.5f %10.5f %10.3f \n" % \
    32133222                                                (int(peak[0]),int(peak[1]),int(peak[2]),int(peak[3]),peak[4],peak[5],FWHM/100.,peak[8]))
  • TabularUnified trunk/GSASIIgrid.py

    r2127 r2129  
    2323import sys
    2424import os
     25import random as ran
    2526import numpy as np
    2627import numpy.ma as ma
     
    29912992            'Scale':1.0,'oldxy':[],'viewDir':[1,0,0]},'Super':Super,'SuperVec':SuperVec}
    29922993        G2plt.Plot3DSngl(G2frame,newPlot=True,Data=controls,hklRef=refList,Title=phaseName)
    2993        
    2994 #    def OnMerge(self,event):
    2995 #        if not len(self.HKL):
    2996 #            print 'No data'
    2997 #            return
    2998 #        self.newHKL = np.copy(self.HKL)
    2999 #        for H in self.newHKL:
    3000 #            H[:4] = np.rint(np.inner(self.Trans,H[:4]))
    3001 #        self.newHKL = np.asarray(self.newHKL)
    3002 #        self.newHKL = G2lat.LaueUnique(self.Laue,self.newHKL)
    3003 #        dlg = wx.ProgressDialog('Build HKL dictonary','',len(self.newHKL)+1,
    3004 #            style = wx.PD_ELAPSED_TIME|wx.PD_AUTO_HIDE)
    3005 #        HKLdict = {}
    3006 #        for ih,hkl in enumerate(self.newHKL):
    3007 #            if str(hkl[:4]) not in HKLdict:
    3008 #                HKLdict[str(hkl[:4])] = [hkl[:4],[hkl[4:],]]
    3009 #            else:
    3010 #                HKLdict[str(hkl[:4])][1].append(hkl[4:])
    3011 #            dlg.Update(ih)
    3012 #        dlg.Destroy()
    3013 #        self.newHKL = []
    3014 #        dlg = wx.ProgressDialog('Processing merge','',len(HKLdict)+1,
    3015 #            style = wx.PD_ELAPSED_TIME|wx.PD_AUTO_HIDE)
    3016 #        sumDf = 0.
    3017 #        sumFo = 0.
    3018 #        for ih,hkl in enumerate(HKLdict):
    3019 #            HKL = HKLdict[hkl]
    3020 #            newHKL = list(HKL[0])
    3021 #            if len(HKL[1]) > 1:
    3022 #                fos = np.array(HKL[1])
    3023 #                wFo = 1/fos[:,1]**2
    3024 #                Fo = np.average(fos[:,0],weights=wFo)
    3025 #                std = np.std(fos[:,0])
    3026 #                sig = np.sqrt(np.mean(fos[:,1])**2+std**2)
    3027 #                sumFo += np.sum(fos[:,0])
    3028 #                sumDf += np.sum(np.abs(fos[:,0]-Fo))
    3029 #                dlg.Update(ih)
    3030 #            else:
    3031 #                Fo = HKL[1][0][0]
    3032 #                sig = HKL[1][0][1]
    3033 #            newHKL += [Fo,sig]
    3034 #            if Fo > 0.:
    3035 #                self.newHKL.append(list(newHKL))
    3036 #        dlg.Destroy()
    3037 #        self.newHKL = np.array(self.newHKL)
    3038 #        print 'merge R = %6.2f%s for %d reflections'%(100.*sumDf/sumFo,'%',self.newHKL.shape[0])
    3039 #        self.newHKL = G2mth.sortArray(G2mth.sortArray(G2mth.sortArray(G2mth.sortArray(self.newHKL,3),2),1),0)
    3040            
     2994                 
    30412995    def OnMergeHKL(event):
     2996        Name = G2frame.PatternTree.GetItemText(G2frame.PatternId)
     2997        Inst = G2frame.PatternTree.GetItemPyData(GetPatternTreeItemId(G2frame,
     2998            G2frame.PatternId,'Instrument Parameters'))
     2999        refList = np.copy(data[1]['RefList'])
    30423000        dlg = MergeDialog(G2frame,data)
    30433001        try:
    30443002            if dlg.ShowModal() == wx.ID_OK:
    30453003                Trans,Cent,Laue = dlg.GetSelection()
    3046                 print "do merge here: ",Cent,Laue
    30473004            else:
    30483005                return
    30493006        finally:
    30503007            dlg.Destroy()
    3051        
    3052         print 'merge HKLF'
    3053            
     3008        Super = data[1]['Super']
     3009        refList = G2lat.transposeHKLF(Trans,Super,refList)
     3010        refList = G2lat.LaueUnique(Laue,refList)
     3011        dlg = wx.ProgressDialog('Build HKL dictonary','',len(refList)+1,
     3012            style = wx.PD_ELAPSED_TIME|wx.PD_AUTO_HIDE)
     3013        HKLdict = {}
     3014        for ih,hkl in enumerate(refList):
     3015            if str(hkl[:3+Super]) not in HKLdict:
     3016                HKLdict[str(hkl[:3+Super])] = [hkl[:3+Super],[hkl[3+Super:],]]
     3017            else:
     3018                HKLdict[str(hkl[:3+Super])][1].append(hkl[3+Super:])
     3019            dlg.Update(ih)
     3020        dlg.Destroy()
     3021        mergeRef = []
     3022        dlg = wx.ProgressDialog('Processing merge','',len(HKLdict)+1,
     3023            style = wx.PD_ELAPSED_TIME|wx.PD_AUTO_HIDE)
     3024        sumDf = 0.
     3025        sumFo = 0.
     3026        for ih,hkl in enumerate(HKLdict):
     3027            HKL = HKLdict[hkl]
     3028            newHKL = list(HKL[0])+list(HKL[1][0])
     3029            if len(HKL[1]) > 1:
     3030                fos = np.array(HKL[1])
     3031                wFo = 1/fos[:,3]**2
     3032                Fo = np.average(fos[:,2],weights=wFo)
     3033                std = np.std(fos[:,2])
     3034                sig = np.sqrt(np.mean(fos[:,3])**2+std**2)
     3035                sumFo += np.sum(fos[:,2])
     3036                sumDf += np.sum(np.abs(fos[:,2]-Fo))
     3037                dlg.Update(ih)
     3038                newHKL[5+Super] = Fo
     3039                newHKL[6+Super] = sig
     3040                newHKL[8+Super] = Fo
     3041            if newHKL[5+Super] > 0.:
     3042                mergeRef.append(list(newHKL))
     3043        dlg.Destroy()
     3044        if Super:
     3045            mergeRef = G2mth.sortArray(G2mth.sortArray(G2mth.sortArray(G2mth.sortArray(mergeRef,3),2),1),0)
     3046        else:
     3047            mergeRef = G2mth.sortArray(G2mth.sortArray(G2mth.sortArray(mergeRef,2),1),0)
     3048        mergeRef = np.array(mergeRef)
     3049        print 'merge R = %6.2f%s for %d reflections'%(100.*sumDf/sumFo,'%',mergeRef.shape[0])
     3050        HKLFlist = []
     3051        newName = Name+' '+Laue
     3052        if G2frame.PatternTree.GetCount():
     3053            item, cookie = G2frame.PatternTree.GetFirstChild(G2frame.root)
     3054            while item:
     3055                name = G2frame.PatternTree.GetItemText(item)
     3056                if name.startswith('HKLF ') and name not in HKLFlist:
     3057                    HKLFlist.append(name)
     3058                item, cookie = G2frame.PatternTree.GetNextChild(G2frame.root, cookie)
     3059        newName = G2obj.MakeUniqueLabel(newName,HKLFlist)
     3060        newData = copy.deepcopy(data)
     3061        newData[0]['ranId'] = ran.randint(0,sys.maxint)
     3062        newData[1]['RefList'] = mergeRef
     3063        Id = G2frame.PatternTree.AppendItem(parent=G2frame.root,text=newName)
     3064        G2frame.PatternTree.SetItemPyData(Id,newData)
     3065        G2frame.PatternTree.SetItemPyData(
     3066            G2frame.PatternTree.AppendItem(Id,text='Instrument Parameters'),Inst)
     3067        G2frame.PatternTree.SetItemPyData(
     3068            G2frame.PatternTree.AppendItem(Id,text='Reflection List'),{})  #dummy entry for GUI use
     3069                   
    30543070    def OnErrorAnalysis(event):
    30553071        G2plt.PlotDeltSig(G2frame,kind)
  • TabularUnified trunk/GSASIIlattice.py

    r2126 r2129  
    524524    return Hmax
    525525   
     526def transposeHKLF(transMat,Super,refList):
     527    newRefs = np.copy(refList)
     528    for H in newRefs:
     529        H[:3+Super] = np.rint(np.inner(transMat,H[:3+Super]))
     530    return newRefs
     531   
    526532def sortHKLd(HKLd,ifreverse,ifdup,ifSS=False):
    527     '''needs doc string
     533    '''sort reflection list on d-spacing; can sort in either order
    528534
    529535    :param HKLd: a list of [h,k,l,d,...];
     
    940946    mat4 = np.array([[0,1,0],[-1,0,0],[0,0,1]])
    941947    mat3 = np.array([[0,-1,0],[1,-1,0],[0,0,1]])
     948    mat6 = np.array([[0,1,0],[-1,1,0],[0,0,1]])
    942949    #triclinic
    943     if Laue == '1':
     950    if Laue == '1': #ok
    944951        pass
    945     elif Laue == '-1':
     952    elif Laue == '-1':  #ok
    946953        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,-1])[:,nxs],HKLFT[:3])
     954        HKLFT[:3] = np.where(HKLFT[0]==0,np.where(HKLFT[1]<0,HKLFT[:3]*np.array([-1,-1,-1])[:,nxs],HKLFT[:3]),HKLFT[:3])
    947955    #monoclinic - all 3 settings
    948     elif Laue == '2(a)':       
    949         HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([1,-1,-1])[:,nxs],HKLFT[:3])
    950     elif Laue == '2':
    951         HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,1,-1])[:,nxs],HKLFT[:3])
    952     elif Laue == '2(c)':       
    953         HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,1])[:,nxs],HKLFT[:3])
    954     elif Laue == 'm(a)':       
     956    elif Laue == '2(a)':    #ok 
     957        HKLFT[:3] = np.where(HKLFT[2]<0,HKLFT[:3]*np.array([1,-1,-1])[:,nxs],HKLFT[:3])
     958        HKLFT[:3] = np.where(HKLFT[2]==0,np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,1])[:,nxs],HKLFT[:3]),HKLFT[:3])
     959    elif Laue == '2':   #ok
     960        HKLFT[:3] = np.where(HKLFT[0]<=0,HKLFT[:3]*np.array([-1,1,-1])[:,nxs],HKLFT[:3])
     961        HKLFT[:3] = np.where(HKLFT[0]==0,np.where(HKLFT[1]<0,HKLFT[:3]*np.array([1,-1,-1])[:,nxs],HKLFT[:3]),HKLFT[:3])
     962    elif Laue == '2(c)':   #ok   
     963        HKLFT[:3] = np.where(HKLFT[1]<0,HKLFT[:3]*np.array([-1,-1,1])[:,nxs],HKLFT[:3])
     964        HKLFT[:3] = np.where(HKLFT[1]==0,np.where(HKLFT[2]<0,HKLFT[:3]*np.array([-1,1,-1])[:,nxs],HKLFT[:3]),HKLFT[:3])
     965    elif Laue == 'm(a)':        #ok   
    955966        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,1,1])[:,nxs],HKLFT[:3])
    956     elif Laue == 'm':
    957         HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([1,-1,1])[:,nxs],HKLFT[:3])
    958     elif Laue == 'm(c)':       
    959         HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([1,1,-1])[:,nxs],HKLFT[:3])
     967    elif Laue == 'm':           #ok
     968        HKLFT[:3] = np.where(HKLFT[1]<0,HKLFT[:3]*np.array([1,-1,1])[:,nxs],HKLFT[:3])
     969    elif Laue == 'm(c)':        #ok
     970        HKLFT[:3] = np.where(HKLFT[2]<0,HKLFT[:3]*np.array([1,1,-1])[:,nxs],HKLFT[:3])
    960971    elif Laue == '2/m(a)':       
    961972        HKLFT[:3] = np.where(HKLFT[1]<0,HKLFT[:3]*np.array([1,-1,-1])[:,nxs],HKLFT[:3])
    962         HKLFT[:3] = np.where(HKLFT[2]<0,HKLFT[:3]*np.array([-1,1,1])[:,nxs],HKLFT[:3])
    963     elif Laue == '2/m':
     973        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,1,1])[:,nxs],HKLFT[:3])
     974    elif Laue == '2/m': #ok
    964975        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,1,-1])[:,nxs],HKLFT[:3])
    965         HKLFT[:3] = np.where(HKLFT[2]<0,HKLFT[:3]*np.array([1,-1,1])[:,nxs],HKLFT[:3])
     976        HKLFT[:3] = np.where(HKLFT[1]<0,HKLFT[:3]*np.array([1,-1,1])[:,nxs],HKLFT[:3])
    966977    elif Laue == '2/m(c)':
    967         HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,1])[:,nxs],HKLFT[:3])
    968         HKLFT[:3] = np.where(HKLFT[1]<0,HKLFT[:3]*np.array([1,1,-1])[:,nxs],HKLFT[:3])
     978        HKLFT[:3] = np.where(HKLFT[1]<0,HKLFT[:3]*np.array([-1,-1,1])[:,nxs],HKLFT[:3])
     979        HKLFT[:3] = np.where(HKLFT[2]<0,HKLFT[:3]*np.array([1,1,-1])[:,nxs],HKLFT[:3])
    969980    #orthorhombic - 3 settings
    970     elif Laue == '222':
     981    elif Laue == '222': #ok
    971982        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,1])[:,nxs],HKLFT[:3])
    972983        HKLFT[:3] = np.where(HKLFT[1]<0,HKLFT[:3]*np.array([1,-1,-1])[:,nxs],HKLFT[:3])
    973     elif Laue == '2mm':       
    974         HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([1,1,-1])[:,nxs],HKLFT[:3])
    975         HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([1,-1,1])[:,nxs],HKLFT[:3])
    976     elif Laue == 'm2m':       
     984        HKLFT[:3] = np.where(HKLFT[0]==0,np.where(HKLFT[2]<0,HKLFT[:3]*np.array([1,1,-1])[:,nxs],HKLFT[:3]),HKLFT[:3])
     985        HKLFT[:3] = np.where(HKLFT[1]==0,np.where(HKLFT[2]<0,HKLFT[:3]*np.array([1,1,-1])[:,nxs],HKLFT[:3]),HKLFT[:3])
     986    elif Laue == '2mm':    #ok   
     987        HKLFT[:3] = np.where(HKLFT[1]<0,HKLFT[:3]*np.array([1,-1,1])[:,nxs],HKLFT[:3])
     988        HKLFT[:3] = np.where(HKLFT[2]<0,HKLFT[:3]*np.array([1,1,-1])[:,nxs],HKLFT[:3])
     989    elif Laue == 'm2m':       #ok
    977990        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,1,1])[:,nxs],HKLFT[:3])
    978         HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([1,1,-1])[:,nxs],HKLFT[:3])
    979     elif Laue == 'mm2':       
     991        HKLFT[:3] = np.where(HKLFT[2]<0,HKLFT[:3]*np.array([1,1,-1])[:,nxs],HKLFT[:3])
     992    elif Laue == 'mm2':    #ok   
    980993        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,1,1])[:,nxs],HKLFT[:3])
    981994        HKLFT[:3] = np.where(HKLFT[1]<0,HKLFT[:3]*np.array([1,-1,1])[:,nxs],HKLFT[:3])
    982     elif Laue == 'mmm':
     995    elif Laue == 'mmm': #ok
    983996        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,1,1])[:,nxs],HKLFT[:3])
    984997        HKLFT[:3] = np.where(HKLFT[1]<0,HKLFT[:3]*np.array([1,-1,1])[:,nxs],HKLFT[:3])
     
    9871000    elif Laue == '-3':
    9881001        HKLFT[:3] = np.where(HKLFT[2]<0,HKLFT[:3]*np.array([-1,-1,-1])[:,nxs],HKLFT[:3])
    989         HKLFT[:3] = np.where(HKLFT[0]<HKLF[1],np.inner(HKLFT[:3]*mat3[:,:,nxs]),HKLFT[:3])
     1002        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat3[nxs,:,:])).T,HKLFT[:3])
    9901003    elif Laue == '3/m':
    9911004        HKLFT[:3] = np.where(HKLFT[2]<0,HKLFT[:3]*np.array([1,1,-1])[:,nxs],HKLFT[:3])
    992         HKLFT[:3] = np.where(HKLFT[0]<HKLF[1],np.inner(HKLFT[:3]*mat3[:,:,nxs]),HKLFT[:3])
     1005        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat3[nxs,:,:])).T,HKLFT[:3])
    9931006    elif Laue == '3':
    994         HKLFT[:3] = np.where(HKLFT[0]<HKLF[1],np.inner(HKLFT[:3]*mat3[:,:,nxs]),HKLFT[:3])
     1007        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat3[nxs,:,:])).T,HKLFT[:3])
    9951008    elif Laue == '32':
    9961009        HKLFT[:3] = np.where(HKLFT[2]<0,HKLFT[:3]*np.array([-1,1,-1])[:,nxs],HKLFT[:3])
    997         HKLFT[:3] = np.where(HKLFT[0]<HKLF[1],np.inner(HKLFT[:3]*mat3[:,:,nxs]),HKLFT[:3])
    998         HKLFT[:3] = np.where(HKLFT[0]<HKLF[1],np.inner(HKLFT[:3]*mat3[:,:,nxs]),HKLFT[:3])
     1010        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat3[nxs,:,:])).T,HKLFT[:3])
    9991011    elif Laue == '3m':
    10001012        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,1,1])[:,nxs],HKLFT[:3])
    1001         HKLFT[:3] = np.where(HKLFT[0]<HKLF[1],np.inner(HKLFT[:3]*mat3[:,:,nxs]),HKLFT[:3])
     1013        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat3[nxs,:,:])).T,HKLFT[:3])
    10021014    #tetragonal
    10031015    elif Laue == '4/m':
    10041016        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,1])[:,nxs],HKLFT[:3])
    10051017        HKLFT[:3] = np.where(HKLFT[1]<0,HKLFT[:3]*np.array([1,1,-1])[:,nxs],HKLFT[:3])
    1006         HKLFT[:3] = np.where(HKLFT[1]<HKLF[0],np.inner(HKLFT[:3]*mat4[:,:,nxs]),HKLFT[:3])
     1018        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat4[nxs,:,:])).T,HKLFT[:3])
    10071019    elif Laue == '4/mmm':
    10081020        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,-1])[:,nxs],HKLFT[:3])
    10091021        HKLFT[:3] = np.where(HKLFT[1]<0,HKLFT[:3]*np.array([-1,-1,-1])[:,nxs],HKLFT[:3])
    10101022        HKLFT[:3] = np.where(HKLFT[2]<0,HKLFT[:3]*np.array([-1,-1,-1])[:,nxs],HKLFT[:3])
    1011         HKLFT[:3] = np.where(HKLFT[1]<HKLF[0],np.inner(HKLFT[:3]*mat4[:,:,nxs]),HKLFT[:3])
     1023        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat4[nxs,:,:])).T,HKLFT[:3])
    10121024    elif Laue == '4':
    1013         HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,1])[:,nxs],HKLFT[:3])
    1014         HKLFT[:3] = np.where(HKLFT[1]<HKLF[0],np.inner(HKLFT[:3]*mat4[:,:,nxs]),HKLFT[:3])
     1025        HKLFT[:3] = np.where(HKLFT[1]<0,HKLFT[:3]*np.array([-1,-1,1])[:,nxs],HKLFT[:3])
     1026        HKLFT[:3] = np.where(HKLFT[1]==0,np.where(HKLFT[2]<0,HKLFT[:3]*np.array([-1,1,-1])[:,nxs],HKLFT[:3]),HKLFT[:3])
     1027        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat4[nxs,:,:])).T,HKLFT[:3])
    10151028    elif Laue == '-4':
    10161029        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,-1])[:,nxs],HKLFT[:3])
    1017         HKLFT[:3] = np.where(HKLFT[1]<HKLF[0],np.inner(HKLFT[:3]*mat4[:,:,nxs]),HKLFT[:3])
     1030        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat4[nxs,:,:])).T,HKLFT[:3])
    10181031    elif Laue == '422':
    10191032        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,1])[:,nxs],HKLFT[:3])
    10201033        HKLFT[:3] = np.where(HKLFT[1]<0,HKLFT[:3]*np.array([1,-1,-1])[:,nxs],HKLFT[:3])
    1021         HKLFT[:3] = np.where(HKLFT[1]<HKLF[0],np.inner(HKLFT[:3]*mat4[:,:,nxs]),HKLFT[:3])
     1034        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat4[nxs,:,:])).T,HKLFT[:3])
    10221035    elif Laue == '-42m':
    10231036        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,-1])[:,nxs],HKLFT[:3])
     1037        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat4[nxs,:,:])).T,HKLFT[:3])
    10241038    elif Laue == '42m':
    10251039        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,-1])[:,nxs],HKLFT[:3])
     1040        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat4[nxs,:,:])).T,HKLFT[:3])
    10261041    #hexagonal
    10271042    elif Laue == '6/m':
    10281043        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,-1])[:,nxs],HKLFT[:3])
     1044        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat3[nxs,:,:])).T,HKLFT[:3])
    10291045    elif Laue == '6/mmm':
    10301046        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,-1])[:,nxs],HKLFT[:3])
     1047        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat3[nxs,:,:])).T,HKLFT[:3])
    10311048    elif Laue == '6':
    10321049        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,-1])[:,nxs],HKLFT[:3])
     1050        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat3[nxs,:,:])).T,HKLFT[:3])
    10331051    elif Laue == '-6':
    10341052        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,-1])[:,nxs],HKLFT[:3])
     1053        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat3[nxs,:,:])).T,HKLFT[:3])
    10351054    elif Laue == '622':
    10361055        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,-1])[:,nxs],HKLFT[:3])
     1056        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat3[nxs,:,:])).T,HKLFT[:3])
    10371057    elif Laue == '-62m':
    10381058        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,-1])[:,nxs],HKLFT[:3])
     1059        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat3[nxs,:,:])).T,HKLFT[:3])
    10391060    elif Laue == '62m':
    10401061        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,-1])[:,nxs],HKLFT[:3])
     1062        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat3[nxs,:,:])).T,HKLFT[:3])
    10411063    #cubic
    10421064    elif Laue == 'm3':
    10431065        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,-1])[:,nxs],HKLFT[:3])
     1066        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat4[nxs,:,:])).T,HKLFT[:3])
    10441067    elif Laue == 'm3m':
    10451068        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,-1])[:,nxs],HKLFT[:3])
     1069        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat4[nxs,:,:])).T,HKLFT[:3])
    10461070    elif Laue == '23':
    10471071        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,-1])[:,nxs],HKLFT[:3])
     1072        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat4[nxs,:,:])).T,HKLFT[:3])
    10481073    elif Laue == '432':
    10491074        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,-1])[:,nxs],HKLFT[:3])
     1075        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat4[nxs,:,:])).T,HKLFT[:3])
    10501076    elif Laue == '-432':
    10511077        HKLFT[:3] = np.where(HKLFT[0]<0,HKLFT[:3]*np.array([-1,-1,-1])[:,nxs],HKLFT[:3])
     1078        HKLFT[:3] = np.where(HKLFT[1]<HKLFT[0],np.squeeze(np.inner(HKLF[:,:3],mat4[nxs,:,:])).T,HKLFT[:3])
    10521079    return HKLFT.T
    10531080       
  • TabularUnified trunk/GSASIImath.py

    r2115 r2129  
    31763176
    31773177def getCWsig(ins,pos):
    3178     '''get CW peak profile sigma
     3178    '''get CW peak profile sigma^2
    31793179   
    31803180    :param dict ins: instrument parameters with at least 'U', 'V', & 'W'
    31813181      as values only
    31823182    :param float pos: 2-theta of peak
    3183     :returns: float getCWsig: peak sigma
     3183    :returns: float getCWsig: peak sigma^2
    31843184   
    31853185    '''
     
    31883188   
    31893189def getCWsigDeriv(pos):
    3190     '''get derivatives of CW peak profile sigma wrt U,V, & W
     3190    '''get derivatives of CW peak profile sigma^2 wrt U,V, & W
    31913191   
    31923192    :param float pos: 2-theta of peak
    31933193   
    3194     :returns: list getCWsigDeriv: d(sig)/dU, d(sig)/dV & d(sig)/dW
     3194    :returns: list getCWsigDeriv: d(sig^2)/dU, d(sig)/dV & d(sig)/dW
    31953195   
    31963196    '''
     
    32203220   
    32213221def getTOFsig(ins,dsp):
    3222     '''get TOF peak profile sigma
     3222    '''get TOF peak profile sigma^2
    32233223   
    32243224    :param dict ins: instrument parameters with at least 'sig-0', 'sig-1' & 'sig-q'
     
    32263226    :param float dsp: d-spacing of peak
    32273227   
    3228     :returns: float getTOFsig: peak sigma
     3228    :returns: float getTOFsig: peak sigma^2
    32293229   
    32303230    '''
     
    32323232   
    32333233def getTOFsigDeriv(dsp):
    3234     '''get derivatives of TOF peak profile gamma wrt sig-0, sig-1, & sig-q
     3234    '''get derivatives of TOF peak profile sigma^2 wrt sig-0, sig-1, & sig-q
    32353235   
    32363236    :param float dsp: d-spacing of peak
     
    40824082            h,k,l,m,d,pos,sig,gam,f = ref[:9]
    40834083            if d >= MCSA['dmin']:
    4084                 sig = G2pwd.getgamFW(sig,gam)/sq8ln2        #--> sig from FWHM
     4084                sig = G2pwd.getgamFW(gam,sig)/sq8ln2        #--> sig from FWHM
    40854085                SQ = 0.25/d**2
    40864086                allFF.append(allM*[G2el.getFFvalues(FFtables,SQ,True)[i] for i in allT]/np.max(allM))
  • TabularUnified trunk/GSASIIpwd.py

    r2049 r2129  
    457457def getWidthsCW(pos,sig,gam,shl):
    458458    '''Compute the peak widths used for computing the range of a peak
    459     for constant wavelength data.
    460     On low-angle side, 10 FWHM are used, on high-angle side 15 are used
     459    for constant wavelength data. On low-angle side, 10 FWHM are used,
     460    on high-angle side 15 are used, low angle side extended by axial divergence
    461461    (for peaks above 90 deg, these are reversed.)
    462462    '''
    463     widths = [np.sqrt(sig)/100.,gam/200.]
    464     fwhm = 2.355*widths[0]+2.*widths[1]
     463    widths = [np.sqrt(sig)/100.,gam/100.]
     464    fwhm = 2.355*widths[0]+widths[1]
    465465    fmin = 10.*(fwhm+shl*abs(npcosd(pos)))
    466466    fmax = 15.0*fwhm
     
    470470   
    471471def getWidthsTOF(pos,alp,bet,sig,gam):
    472     'needs a doc string'
     472    '''Compute the peak widths used for computing the range of a peak
     473    for constant wavelength data. 10 FWHM are used on both sides each
     474    extended by exponential coeff.
     475    '''
    473476    lnf = 3.3      # =log(0.001/2)
    474477    widths = [np.sqrt(sig),gam]
     
    494497   
    495498def getgamFW(g,s):
    496     'needs a doc string'
     499    '''Compute total FWHM from Thompson, Cox & Hastings (1987), J. Appl. Cryst. 20, 79-83
     500   
     501    :param g: float Lorentzian FWHM
     502    :param s: float Gaussian FWHM
     503   
     504    :returns float: total FWHM of pseudoVoigt
     505    '''
    497506    gamFW = lambda s,g: np.exp(np.log(s**5+2.69269*s**4*g+2.42843*s**3*g**2+4.47163*s**2*g**3+0.07842*s*g**4+g**5)/5.)
    498     return gamFW(g,s)
     507    return gamFW(s,g)
    499508               
    500509def getFCJVoigt(pos,intens,sig,gam,shl,xdata):   
     
    901910                else:
    902911                    sig = G2mth.getCWsig(parmDict,tth)
    903                 sig = max(sig,0.001)          #avoid neg sigma
     912                sig = max(sig,0.001)          #avoid neg sigma^2
    904913                gamName = 'gam'+str(iPeak)
    905914                if gamName in varyList:
Note: See TracChangeset for help on using the changeset viewer.