Changeset 390


Ignore:
Timestamp:
Oct 10, 2011 4:32:58 PM (12 years ago)
Author:
vondreele
Message:

Add tooltip to covariance plot
clean up LS output formats
Add Fobs extraction & R-value calculations

Location:
trunk
Files:
3 edited

Legend:

Unmodified
Added
Removed

  • trunk/GSASIIplot.py

    r388 r390  
    13931393            ypos = int(event.ydata+.5)
    13941394            if -1 < xpos < len(varyList) and -1 < ypos < len(varyList):
    1395                 self.G2plotNB.status.SetFields(['Key: s to change colors',
    1396                     '%s - %s: %5.3f'%(varyList[xpos].ljust(19),varyList[ypos].ljust(19),covArray[xpos][ypos])])
     1395                msg = '%s - %s: %5.3f'%(varyList[xpos],varyList[ypos],covArray[xpos][ypos])
     1396                Page.canvas.SetToolTipString(msg)
     1397                self.G2plotNB.status.SetFields(['Key: s to change colors',msg])
    13971398    try:
    13981399        plotNum = self.G2plotNB.plotList.index('Covariance')

  • trunk/GSASIIpwdGUI.py

    r380 r390  
    12931293    rowLabels = []
    12941294    refList = []
    1295     for h,k,l,m,d,pos,sig,gam,fo,fc,phi,x,x in data[self.RefList]:
     1295    for h,k,l,m,d,pos,sig,gam,fo,fc,phi,x,x,x in data[self.RefList]:
    12961296        refList.append([h,k,l,m,d,pos,sig,gam,fo,fc,phi])
    12971297    for i in range(len(refList)): rowLabels.append(str(i))

  • trunk/GSASIIstruct.py

    r386 r390  
    189189       
    190190    GPXback = GPXBackup(GPXfile)
    191     print '\n',130*'-'
     191    print '\n',135*'-'
    192192    print 'Read from file:',GPXback
    193193    print 'Save to file  :',GPXfile
    … …  
    631631   
    632632    def PrintSize(hapData):
    633         line = '\n Size model    : '+hapData[0]
    634633        if hapData[0] in ['isotropic','uniaxial']:
     634            line = '\n Size model    : %9s'%(hapData[0])
    635635            line += ' equatorial:'+'%12.3f'%(hapData[1][0])+' Refine? '+str(hapData[2][0])
    636636            if hapData[0] == 'uniaxial':
    … …  
    638638            print line
    639639        else:
    640             print line
     640            print '\n Size model    : %s'%(hapData[0])
    641641            Snames = ['S11','S22','S33','S12','S13','S23']
    642642            ptlbls = ' names :'
    … …  
    652652       
    653653    def PrintMuStrain(hapData,SGData):
    654         line = '\n Mustrain model: '+hapData[0]
    655654        if hapData[0] in ['isotropic','uniaxial']:
     655            line = '\n Mustrain model: %9s'%(hapData[0])
    656656            line += ' equatorial:'+'%12.1f'%(hapData[1][0])+' Refine? '+str(hapData[2][0])
    657657            if hapData[0] == 'uniaxial':
    … …  
    659659            print line
    660660        else:
    661             print line
     661            print '\n Mustrain model: %s'%(hapData[0])
    662662            Snames = G2spc.MustrainNames(SGData)
    663663            ptlbls = ' names :'
    … …  
    775775            if Print:
    776776                print '\n Phase: ',phase,' in histogram: ',histogram
    777                 print '\n Phase fraction  : %10.4f'%(hapData['Scale'][0]),' Refine?',hapData['Scale'][1]
     777                print 135*'-'
     778                print ' Phase fraction  : %10.4f'%(hapData['Scale'][0]),' Refine?',hapData['Scale'][1]
    778779                print ' Extinction coeff: %10.4f'%(hapData['Extinction'][0]),' Refine?',hapData['Extinction'][1]
    779780                if hapData['Pref.Ori.'][0] == 'MD':
    … …  
    795796                    pos = 2.0*asind(wave/(2.0*d))
    796797                    if limits[0] < pos < limits[1]:
    797                         refList.append([h,k,l,mul,d,pos,0.0,0.0,0.0,0.0,0.0,Uniq,phi])
     798                        refList.append([h,k,l,mul,d,pos,0.0,0.0,0.0,0.0,0.0,Uniq,phi,0.0])
    798799                else:
    799800                    raise ValueError
    … …  
    804805   
    805806    def PrintSizeAndSig(hapData,sizeSig):
    806         line = '\n Size model:     '+hapData[0]
     807        line = '\n Size model:     %9s'%(hapData[0])
    807808        if hapData[0] in ['isotropic','uniaxial']:
    808809            line += ' equatorial:%12.3f'%(hapData[1][0])
    … …  
    832833       
    833834    def PrintMuStrainAndSig(hapData,mustrainSig,SGData):
    834         line = '\n Mustrain model: '+hapData[0]
     835        line = '\n Mustrain model: %9s'%(hapData[0])
    835836        if hapData[0] in ['isotropic','uniaxial']:
    836837            line += ' equatorial:%12.1f'%(hapData[1][0])
    … …  
    898899        for histogram in HistoPhase:
    899900            print '\n Phase: ',phase,' in histogram: ',histogram
     901            print 130*'-'
    900902            hapData = HistoPhase[histogram]
    901903            Histogram = Histograms[histogram]
    902904            hId = Histogram['hId']
    903905            pfx = str(pId)+':'+str(hId)+':'
     906            print ' Final refinement RF, RF^2 = %.2f%%, %.2f%% on %d reflections'   \
     907                %(Histogram[pfx+'Rf'],Histogram[pfx+'Rf^2'],Histogram[pfx+'Nref'])
    904908           
    905909            PhFrExtPOSig = {}
    … …  
    917921                    if pfx+item in sigDict:
    918922                        PhFrExtPOSig[item] = sigDict[pfx+item]
    919             print '\n'
    920923            if 'Scale' in PhFrExtPOSig:
    921924                print ' Phase fraction  : %10.4f, sig %10.4f'%(hapData['Scale'][0],PhFrExtPOSig['Scale'])
    … …  
    12051208            print 135*'-'
    12061209            print ' Instrument type: ',Sample['Type']
     1210            print ' Final refinement wRp = %.2f%% on %d observations in this histogram'%(Histogram['wRp'],Histogram['Nobs'])
    12071211            PrintSampleParmsSig(Sample,sampSig)
    12081212            PrintInstParmsSig(Inst,instSig)
    … …  
    14581462    Icorr = parmDict[phfx+'Scale']*parmDict[hfx+'Scale']*refl[3]               #scale*multiplicity
    14591463    Icorr *= G2pwd.Polarization(parmDict[hfx+'Polariz.'],refl[5],parmDict[hfx+'Azimuth'])[0]
    1460     Icorr *= GetPrefOri(refl,G,g,phfx,hfx,SGData,calcControls,parmDict)       
    1461     return Icorr
     1464    Icorr *= GetPrefOri(refl,G,g,phfx,hfx,SGData,calcControls,parmDict)
     1465    refl[13] = Icorr       
    14621466   
    14631467def GetIntensityDerv(refl,G,g,phfx,hfx,SGData,calcControls,parmDict):
    1464     Icorr = GetIntensityCorr(refl,G,g,phfx,hfx,SGData,calcControls,parmDict)
    14651468    pola,dIdPola = G2pwd.Polarization(parmDict[hfx+'Polariz.'],refl[5],parmDict[hfx+'Azimuth'])
    14661469    POcorr,dIdPO = GetPrefOriDerv(refl,G,g,phfx,hfx,SGData,calcControls,parmDict)
    … …  
    14701473    dIdsh = 1./parmDict[hfx+'Scale']
    14711474    dIdsp = 1./parmDict[phfx+'Scale']
    1472     return Icorr,dIdsh,dIdsp,dIdPola,dIdPO
     1475    return dIdsh,dIdsp,dIdPola,dIdPO
    14731476       
    14741477def GetSampleGam(refl,wave,G,phfx,calcControls,parmDict,sizeEllipse):
    … …  
    17111714                refl[5] += GetHStrainShift(refl,SGData,phfx,parmDict)               #apply hydrostatic strain shift
    17121715                refl[6:8] = GetReflSIgGam(refl,wave,G,hfx,phfx,calcControls,parmDict,sizeEllipse)    #peak sig & gam
    1713                 Icorr = GetIntensityCorr(refl,G,g,phfx,hfx,SGData,calcControls,parmDict)*Vst*Lorenz
     1716                GetIntensityCorr(refl,G,g,phfx,hfx,SGData,calcControls,parmDict)    #puts corrections in refl[13]
     1717                refl[13] *= Vst*Lorenz
    17141718                if 'Pawley' in Phase['General']['Type']:
    17151719                    try:
    … …  
    17251729                elif not iBeg-iFin:     #peak above high limit - done
    17261730                    break
    1727                 yc[iBeg:iFin] += Icorr*refl[9]*G2pwd.getFCJVoigt3(refl[5],refl[6],refl[7],shl,x[iBeg:iFin])    #>90% of time spent here
     1731                yc[iBeg:iFin] += refl[13]*refl[9]*G2pwd.getFCJVoigt3(refl[5],refl[6],refl[7],shl,x[iBeg:iFin])    #>90% of time spent here
    17281732                if Ka2:
    17291733                    pos2 = refl[5]+lamRatio*tand(refl[5]/2.0)       # + 360/pi * Dlam/lam * tan(th)
    … …  
    17351739                    elif not iBeg-iFin:     #peak above high limit - done
    17361740                        return yc,yb
    1737                     yc[iBeg:iFin] += Icorr*refl[9]*kRatio*G2pwd.getFCJVoigt3(pos2,refl[6],refl[7],shl,x[iBeg:iFin])        #and here
     1741                    yc[iBeg:iFin] += refl[13]*refl[9]*kRatio*G2pwd.getFCJVoigt3(pos2,refl[6],refl[7],shl,x[iBeg:iFin])        #and here
    17381742            elif 'T' in calcControls[hfx+'histType']:
    17391743                print 'TOF Undefined at present'
    1740                 raise ValueError    #no TOF yet
    1741     return yc,yb   
    1742            
     1744                raise Exception    #no TOF yet
     1745    return yc,yb
     1746   
     1747def GetFobsSq(Histograms,Phases,parmDict,calcControls):
     1748    for histogram in Histograms:
     1749        if 'PWDR' in histogram[:4]:
     1750            Histogram = Histograms[histogram]
     1751            hId = Histogram['hId']
     1752            hfx = ':%d:'%(hId)
     1753            Limits = calcControls[hfx+'Limits']
     1754            shl = max(parmDict[hfx+'SH/L'],0.002)
     1755            Ka2 = False
     1756            kRatio = 0.0
     1757            if hfx+'Lam1' in parmDict.keys():
     1758                Ka2 = True
     1759                lamRatio = 360*(parmDict[hfx+'Lam2']-parmDict[hfx+'Lam1'])/(np.pi*parmDict[hfx+'Lam1'])
     1760                kRatio = parmDict[hfx+'I(L2)/I(L1)']
     1761            x,y,w,yc,yb,yd = Histogram['Data']
     1762            ymb = np.array(y-yb)
     1763            ycmb = np.array(yc-yb)
     1764            ratio = np.where(ycmb>0.,ymb/ycmb,1.0)           
     1765            xB = np.searchsorted(x,Limits[0])
     1766            xF = np.searchsorted(x,Limits[1])
     1767            refLists = Histogram['Reflection Lists']
     1768            for phase in refLists:
     1769                Phase = Phases[phase]
     1770                pId = Phase['pId']
     1771                phfx = '%d:%d:'%(pId,hId)
     1772                refList = refLists[phase]
     1773                sumFo = 0.0
     1774                sumdF = 0.0
     1775                sumFosq = 0.0
     1776                sumdFsq = 0.0
     1777                for refl in refList:
     1778                    if 'C' in calcControls[hfx+'histType']:
     1779                        yp = np.zeros_like(yb)
     1780                        Wd,fmin,fmax = G2pwd.getWidths(refl[5],refl[6],refl[7],shl)
     1781                        iBeg = np.searchsorted(x[xB:xF],refl[5]-fmin)
     1782                        iFin = np.searchsorted(x[xB:xF],refl[5]+fmax)
     1783                        iFin2 = iFin
     1784                        yp[iBeg:iFin] = refl[13]*refl[9]*G2pwd.getFCJVoigt3(refl[5],refl[6],refl[7],shl,x[iBeg:iFin])    #>90% of time spent here                           
     1785                        if Ka2:
     1786                            pos2 = refl[5]+lamRatio*tand(refl[5]/2.0)       # + 360/pi * Dlam/lam * tan(th)
     1787                            Wd,fmin,fmax = G2pwd.getWidths(pos2,refl[6],refl[7],shl)
     1788                            iBeg2 = np.searchsorted(x,pos2-fmin)
     1789                            iFin2 = np.searchsorted(x,pos2+fmax)
     1790                            yp[iBeg2:iFin2] += refl[13]*refl[9]*kRatio*G2pwd.getFCJVoigt3(pos2,refl[6],refl[7],shl,x[iBeg2:iFin2])        #and here
     1791                        refl[8] = np.sum(yp[iBeg:iFin2]*ratio[iBeg:iFin2])/(refl[13]*(1.+kRatio))
     1792                    elif 'T' in calcControls[hfx+'histType']:
     1793                        print 'TOF Undefined at present'
     1794                        raise Exception    #no TOF yet
     1795                    Fo = np.sqrt(np.abs(refl[8]))
     1796                    Fc = np.sqrt(np.abs(refl[9]))
     1797                    sumFo += Fo
     1798                    sumFosq += refl[8]**2
     1799                    sumdF += np.abs(Fo-Fc)
     1800                    sumdFsq += (refl[8]-refl[9])**2
     1801                Histogram[phfx+'Rf'] = min(100.,(sumdF/sumFo)*100.)
     1802                Histogram[phfx+'Rf^2'] = min(100.,np.sqrt(sumdFsq/sumFosq)*100.)
     1803                Histogram[phfx+'Nref'] = len(refList)
     1804               
    17431805def getPowderProfileDerv(parmDict,x,varylist,Histogram,Phases,calcControls,pawleyLookup):
    17441806   
    … …  
    17981860        A = [parmDict[pfx+'A%d'%(i)] for i in range(6)]
    17991861        G,g = G2lat.A2Gmat(A)       #recip & real metric tensors
    1800         Vst = np.sqrt(nl.det(G))    #V*
    18011862        if 'Pawley' not in Phase['General']['Type']:
    18021863            dFdvDict = StructureFactorDerv(refList,G,hfx,pfx,SGData,calcControls,parmDict)
    … …  
    18071868            if 'C' in calcControls[hfx+'histType']:        #CW powder
    18081869                h,k,l = refl[:3]
    1809                 Lorenz = 1./(2.*sind(refl[5]/2.)**2*cosd(refl[5]/2.))           #Lorentz correction
    1810                 Icorr,dIdsh,dIdsp,dIdpola,dIdPO = GetIntensityDerv(refl,G,g,phfx,hfx,SGData,calcControls,parmDict)
    1811                 Icorr *= Vst*Lorenz
     1870                dIdsh,dIdsp,dIdpola,dIdPO = GetIntensityDerv(refl,G,g,phfx,hfx,SGData,calcControls,parmDict)
    18121871                if 'Pawley' in Phase['General']['Type']:
    18131872                    try:
    … …  
    18291888                dMdipk = G2pwd.getdFCJVoigt3(refl[5],refl[6],refl[7],shl,x[iBeg:iFin])
    18301889                for i in range(1,5):
    1831                     dMdpk[i][iBeg:iFin] += 100.*dx*Icorr*refl[9]*dMdipk[i]
    1832                 dMdpk[0][iBeg:iFin] += 100.*dx*Icorr*refl[9]*dMdipk[0]
     1890                    dMdpk[i][iBeg:iFin] += 100.*dx*refl[13]*refl[9]*dMdipk[i]
     1891                dMdpk[0][iBeg:iFin] += 100.*dx*refl[13]*refl[9]*dMdipk[0]
    18331892                dervDict = {'int':dMdpk[0],'pos':dMdpk[1],'sig':dMdpk[2],'gam':dMdpk[3],'shl':dMdpk[4]}
    18341893                if Ka2:
    … …  
    18401899                        dMdipk2 = G2pwd.getdFCJVoigt3(pos2,refl[6],refl[7],shl,x[iBeg2:iFin2])
    18411900                        for i in range(1,5):
    1842                             dMdpk[i][iBeg2:iFin2] += 100.*dx*Icorr*refl[9]*kRatio*dMdipk2[i]
    1843                         dMdpk[0][iBeg2:iFin2] += 100.*dx*Icorr*refl[9]*kRatio*dMdipk2[0]
    1844                         dMdpk[5][iBeg2:iFin2] += 100.*dx*Icorr*dMdipk2[0]
     1901                            dMdpk[i][iBeg2:iFin2] += 100.*dx*refl[13]*refl[9]*kRatio*dMdipk2[i]
     1902                        dMdpk[0][iBeg2:iFin2] += 100.*dx*refl[13]*refl[9]*kRatio*dMdipk2[0]
     1903                        dMdpk[5][iBeg2:iFin2] += 100.*dx*refl[13]*dMdipk2[0]
    18451904                        dervDict = {'int':dMdpk[0],'pos':dMdpk[1],'sig':dMdpk[2],'gam':dMdpk[3],'shl':dMdpk[4],'L1/L2':dMdpk[5]*refl[9]}
    18461905                if 'Pawley' in Phase['General']['Type']:
    … …  
    19381997                xB = np.searchsorted(x,Limits[0])
    19391998                xF = np.searchsorted(x,Limits[1])
    1940                 Nobs += xF-xB
     1999                Histogram['Nobs'] = xF-xB
     2000                Nobs += Histogram['Nobs']
    19412001                Histogram['sumwYo'] = np.sum(w[xB:xF]*y[xB:xF]**2)
    19422002                sumwYo += Histogram['sumwYo']
    19432003                yc[xB:xF],yb[xB:xF] = getPowderProfile(parmdict,x[xB:xF],
    19442004                    varylist,Histogram,Phases,calcControls,pawleyLookup)
    1945                 yc[xB:xF] *= parmDict[hfx+'Scale']
    19462005                yc[xB:xF] += yb[xB:xF]
    19472006                yd[xB:xF] = yc[xB:xF]-y[xB:xF]          #yc-yo then all dydv have no '-' needed
    19482007                Histogram['sumwYd'] = np.sum(np.sqrt(w[xB:xF])*(yd[xB:xF]))
    1949                 M = np.concatenate((M,np.sqrt(w[xB:xF])*(yd[xB:xF])))
     2008                wdy = np.sqrt(w[xB:xF])*(yd[xB:xF])
     2009                Histogram['wRp'] = min(100.,np.sqrt(np.sum(wdy**2)/Histogram['sumwYo'])*100.)
     2010                M = np.concatenate((M,wdy))
    19502011        Histograms['sumwYo'] = sumwYo
    19512012        Histograms['Nobs'] = Nobs
    19522013        Rwp = min(100.,np.sqrt(np.sum(M**2)/sumwYo)*100.)
    19532014        if dlg:
    1954             GoOn = dlg.Update(Rwp,newmsg='%s%8.3f%s'%('Powder profile Rwp =',Rwp,'%'))[0]
     2015            GoOn = dlg.Update(Rwp,newmsg='%s%8.3f%s'%('Powder profile wRp =',Rwp,'%'))[0]
    19552016            if not GoOn:
    19562017                return -M           #abort!!
    … …  
    20182079        print ' Number of function calls:',result[2]['nfev'],' Number of observations: ',Histograms['Nobs'],' Number of parameters: ',len(varyList)
    20192080        print ' Refinement time = %8.3fs, %8.3fs/cycle'%(runtime,runtime/ncyc)
    2020         print ' Rwp = %7.2f%%, chi**2 = %12.6g, reduced chi**2 = %6.2f'%(Rwp,chisq,GOF)
     2081        print ' wRp = %7.2f%%, chi**2 = %12.6g, reduced chi**2 = %6.2f'%(Rwp,chisq,GOF)
    20212082        try:
    20222083            covMatrix = result[1]*GOF
    … …  
    20432104#    print 'indParmList: ',G2mv.indParmList
    20442105#    print 'fixedDict: ',G2mv.fixedDict
     2106    GetFobsSq(Histograms,Phases,parmDict,calcControls)
    20452107    sigDict = dict(zip(varyList,sig))
    20462108    SetPhaseData(parmDict,sigDict,Phases)
Note: See TracChangeset for help on using the changeset viewer.