Changeset 1460


Ignore:
Timestamp:
Aug 12, 2014 4:08:52 PM (9 years ago)
Author:
vondreele
Message:

fix missing radius Factor for MC/SA plots
add Prfo, Abs & xt to reflection list items
fix most TOF derivatives (remove a multiply by 100!)

Location:
trunk
Files:
4 edited

Legend:

Unmodified
Added
Removed

  • trunk/GSASIIplot.py

    r1459 r1460  
    37233723
    37243724    def FindPeaksBonds(XYZ):
    3725         rFact = data['Drawing']['radiusFactor']
     3725        rFact = data['Drawing'].get('radiusFactor',0.85)    #data['Drawing'] could be empty!
    37263726        Bonds = [[] for x in XYZ]
    37273727        for i,xyz in enumerate(XYZ):

  • trunk/GSASIIpwdGUI.py

    r1459 r1460  
    25162516            I100 *= 100.0/Imax
    25172517        if 'C' in Inst['Type'][0]:
    2518             refs = np.vstack((refList.T[:11],I100)).T
     2518            refs = np.vstack((refList.T[:15],I100)).T
    25192519        elif 'T' in Inst['Type'][0]:
    2520             refs = np.vstack((refList.T[:15],I100)).T
     2520            refs = np.vstack((refList.T[:18],I100)).T
    25212521           
    25222522    for i in range(len(refs)): rowLabels.append(str(i))
    … …  
    25252525        [wg.GRID_VALUE_FLOAT+':10,3',]
    25262526    if HKLF:
    2527         colLabels = ['H','K','L','mul','d','Fosq','sig','Fcsq','FoTsq','FcTsq','phase','Ext',]
     2527        colLabels = ['H','K','L','mul','d','Fosq','sig','Fcsq','FoTsq','FcTsq','phase','ExtC',]
    25282528    else:
    25292529        if 'C' in Inst['Type'][0]:
    2530             colLabels = ['H','K','L','mul','d','pos','sig','gam','Fosq','Fcsq','phase','Icorr','I100',]
    2531             Types += [wg.GRID_VALUE_FLOAT+':10,3',]
     2530            colLabels = ['H','K','L','mul','d','pos','sig','gam','Fosq','Fcsq','phase','Icorr','Prfo','Trans','ExtP','I100']
     2531            Types += 4*[wg.GRID_VALUE_FLOAT+':10,3',]
    25322532        elif 'T' in Inst['Type'][0]:
    2533             colLabels = ['H','K','L','mul','d','pos','sig','gam','Fosq','Fcsq','phase','Icorr','alp','bet','wave','I100',]
    2534             Types += 4*[wg.GRID_VALUE_FLOAT+':10,3',]
     2533            colLabels = ['H','K','L','mul','d','pos','sig','gam','Fosq','Fcsq','phase','Icorr','alp','bet','wave','Prfo','Abs','Ext','I100']
     2534            Types += 7*[wg.GRID_VALUE_FLOAT+':10,3',]
    25352535           
    25362536    G2frame.PeakTable = G2gd.Table(refs,rowLabels=rowLabels,colLabels=colLabels,types=Types)

  • trunk/GSASIIstrIO.py

    r1459 r1460  
    18131813                            pos = 2.0*asind(wave/(2.0*d))+Zero
    18141814                            if limits[0] < pos < limits[1]:
    1815                                 refList.append([h,k,l,mul,d, pos,0.0,0.0,0.0,0.0, 0.0,0.0])
     1815                                refList.append([h,k,l,mul,d, pos,0.0,0.0,0.0,0.0, 0.0,0.0,1.0,1.0,1.0])
     1816                                #... sig,gam,fotsq,fctsq, phase,icorr,prfo,abs,ext
    18161817                                Uniq.append(uniq)
    18171818                                Phi.append(phi)
    … …  
    18201821                            if limits[0] < pos < limits[1]:
    18211822                                wave = inst['difC'][1]*d/(252.816*inst['fltPath'][0])
    1822                                 refList.append([h,k,l,mul,d, pos,0.0,0.0,0.0,0.0, 0.0,0.0,0.0,0.0,wave])
     1823                                refList.append([h,k,l,mul,d, pos,0.0,0.0,0.0,0.0, 0.0,0.0,0.0,0.0,wave, 1.0,1.0,1.0])
     1824                                # ... sig,gam,fotsq,fctsq, phase,icorr,alp,bet,wave, prfo,abs,ext
    18231825                                Uniq.append(uniq)
    18241826                                Phi.append(phi)

  • trunk/GSASIIstrMath.py

    r1459 r1460  
    802802    corrections needed for derivatives
    803803    '''
    804     ref[11] = 1.0
     804    extCor = 1.0
    805805    dervCor = 1.0
    806806    dervDict = {}
    … …  
    869869            dervDict[phfx+'Eg'] = -ref[7]*PLZ*PF3*(PSIG/parmDict[phfx+'Eg'])**3*PL**2
    870870               
    871         ref[11] = 1./extCor
    872     return dervCor,dervDict
     871    return 1./extCor,dervCor,dervDict
    873872   
    874873def Dict2Values(parmdict, varylist):
    … …  
    10181017   
    10191018def GetPrefOri(refl,uniq,G,g,phfx,hfx,SGData,calcControls,parmDict):
    1020     'Needs a doc string'
     1019    'March-Dollase preferred orientation correction'
    10211020    POcorr = 1.0
    1022     if calcControls[phfx+'poType'] == 'MD':
    1023         MD = parmDict[phfx+'MD']
    1024         if MD != 1.0:
    1025             MDAxis = calcControls[phfx+'MDAxis']
    1026             sumMD = 0
    1027             for H in uniq:           
    1028                 cosP,sinP = G2lat.CosSinAngle(H,MDAxis,G)
    1029                 A = 1.0/np.sqrt((MD*cosP)**2+sinP**2/MD)
    1030                 sumMD += A**3
    1031             POcorr = sumMD/len(uniq)
    1032     else:   #spherical harmonics
    1033         if calcControls[phfx+'SHord']:
    1034             POcorr = SHPOcal(refl,g,phfx,hfx,SGData,calcControls,parmDict)
     1021    MD = parmDict[phfx+'MD']
     1022    if MD != 1.0:
     1023        MDAxis = calcControls[phfx+'MDAxis']
     1024        sumMD = 0
     1025        for H in uniq:           
     1026            cosP,sinP = G2lat.CosSinAngle(H,MDAxis,G)
     1027            A = 1.0/np.sqrt((MD*cosP)**2+sinP**2/MD)
     1028            sumMD += A**3
     1029        POcorr = sumMD/len(uniq)
    10351030    return POcorr
    10361031   
    … …  
    10741069            return G2pwd.AbsorbDerv('Cylinder',parmDict[hfx+'Absorption'],refl[5],0,0)
    10751070    else:
    1076         return G2pwd.SurfaceRoughDerv(parmDict[hfx+'SurfRoughA'],parmDict[hfx+'SurfRoughB'],refl[5])
     1071        return np.array(G2pwd.SurfaceRoughDerv(parmDict[hfx+'SurfRoughA'],parmDict[hfx+'SurfRoughB'],refl[5]))
     1072       
     1073def GetPwdrExt(refl,pfx,phfx,hfx,calcControls,parmDict):
     1074    'Needs a doc string'
     1075    coef = np.array([-0.5,0.25,-0.10416667,0.036458333,-0.0109375,2.8497409E-3])
     1076    pi2 = np.sqrt(np.pi/2.)
     1077    if 'T' in calcControls[hfx+'histType']:
     1078        sth2 = sind(parmDict[hfx+'2-theta']/2.)**2
     1079        wave = refl[12]
     1080    else:   #'C'W
     1081        sth2 = sind(refl[5]/2.)**2
     1082        wave = parmDict.get(hfx+'Lam',parmDict.get(hfx+'Lam1',1.0))
     1083    c2th = 1.-2.0*sth2
     1084    flv2 = refl[9]*(wave/parmDict[pfx+'Vol'])**2
     1085    if 'X' in calcControls[hfx+'histType']:
     1086        flv2 *= 0.079411*(1.0+c2th**2)/2.0
     1087    xfac = flv2*parmDict[phfx+'Extinction']
     1088    exb = 1.0
     1089    if xfac > -1.:
     1090        exb = 1./(1.+xfac)
     1091    exl = 1.0
     1092    if 0 < xfac <= 1.:
     1093        xn = np.array([xfac**(i+1) for i in range(6)])
     1094        exl = np.sum(xn*coef)
     1095    elif xfac > 1.:
     1096        xfac2 = 1./np.sqrt(xfac)
     1097        exl = pi2*(1.-0.125/xfac)*xfac2
     1098    return exb*sth2+exl*(1.-sth2)
     1099   
     1100def GetPwdrExtDerv(refl,pfx,phfx,hfx,calcControls,parmDict):
     1101    'Needs a doc string'
     1102    coef = np.array([-0.5,0.25,-0.10416667,0.036458333,-0.0109375,2.8497409E-3])
     1103    pi2 = np.sqrt(np.pi/2.)
     1104    if 'T' in calcControls[hfx+'histType']:
     1105        sth2 = sind(parmDict[hfx+'2-theta']/2.)**2
     1106        wave = refl[12]
     1107    else:   #'C'W
     1108        sth2 = sind(refl[5]/2.)**2
     1109        wave = parmDict.get(hfx+'Lam',parmDict.get(hfx+'Lam1',1.0))
     1110    c2th = 1.-2.0*sth2
     1111    flv2 = refl[9]*(wave/parmDict[pfx+'Vol'])**2
     1112    return 0.
     1113#
     1114#      STH2 = STHETA**2
     1115#      C2TH = 1-2.0*STH2
     1116#      FLV2 = FCSQ*(LAM/VOL(IPHAS))**2
     1117#      IF ( HTYPE(2:2).EQ.'X' ) FLV2 = 0.079411*FLV2*(1.0+C2TH**2)/2.0
     1118#      XFAC = FLV2*EXTPOWD(IHST,IPHAS)
     1119#      IF ( XFAC.LE.-1.0 ) THEN
     1120#        EXB = 1.0
     1121#        DBDE = -500.0*FLV2
     1122#      ELSE
     1123#        EXB = 1.0/SQRT(1.0+XFAC)
     1124#        DBDE = -0.5*FLV2*EXB**3
     1125#      END IF
     1126#      IF ( XFAC.LE.0.0 ) THEN
     1127#        EXL = 1.0
     1128#        DLDE = 0.0
     1129#      ELSE IF ( XFAC.LE.1.0 ) THEN
     1130#        EXL = 1.0
     1131#        DLDE = 0.0
     1132#        DO I=1,6
     1133#          XN =XFAC**I
     1134#          EXL = EXL+COEF(I)*XN
     1135#          DLDE = DLDE+FLOAT(I)*FLV2*COEF(I)*XN/XFAC
     1136#        END DO
     1137#      ELSE
     1138#        XFAC2 = 1.0/SQRT(XFAC)
     1139#        EXL = PI2*(1.0-0.125/XFAC)*XFAC2
     1140#        DLDE = 0.5*FLV2*PI2*XFAC2*(-1.0/XFAC+0.375/XFAC**2)
     1141#      END IF
     1142#      EXTCOR = EXB*STH2+EXL*(1.0-STH2)
     1143#      DFDEX = DBDE*STH2+DLDE*(1.0-STH2)
    10771144   
    10781145def GetIntensityCorr(refl,uniq,G,g,pfx,phfx,hfx,SGData,calcControls,parmDict):
    … …  
    10811148    if 'X' in parmDict[hfx+'Type']:
    10821149        Icorr *= G2pwd.Polarization(parmDict[hfx+'Polariz.'],refl[5],parmDict[hfx+'Azimuth'])[0]
    1083     Icorr *= GetPrefOri(refl,uniq,G,g,phfx,hfx,SGData,calcControls,parmDict)
    1084     if pfx+'SHorder' in parmDict:
    1085         Icorr *= SHTXcal(refl,g,pfx,hfx,SGData,calcControls,parmDict)
    1086     Icorr *= GetAbsorb(refl,hfx,calcControls,parmDict)
    1087     return Icorr
     1150    POcorr = 1.0
     1151    if pfx+'SHorder' in parmDict:                 #generalized spherical harmonics texture
     1152        POcorr = SHTXcal(refl,g,pfx,hfx,SGData,calcControls,parmDict)
     1153    elif calcControls[phfx+'poType'] == 'MD':         #March-Dollase
     1154        POcorr = GetPrefOri(refl,uniq,G,g,phfx,hfx,SGData,calcControls,parmDict)
     1155    elif calcControls[phfx+'SHord']:                #cylindrical spherical harmonics
     1156        POcorr = SHPOcal(refl,g,phfx,hfx,SGData,calcControls,parmDict)
     1157    Icorr *= POcorr
     1158    AbsCorr = 1.0
     1159    AbsCorr = GetAbsorb(refl,hfx,calcControls,parmDict)
     1160    Icorr *= AbsCorr
     1161    ExtCorr = GetPwdrExt(refl,pfx,phfx,hfx,calcControls,parmDict)
     1162    Icorr *= ExtCorr
     1163    return Icorr,POcorr,AbsCorr,ExtCorr
    10881164   
    10891165def GetIntensityDerv(refl,uniq,G,g,pfx,phfx,hfx,SGData,calcControls,parmDict):
    … …  
    10961172    else:       #'N'
    10971173        dIdPola = 0.0
    1098     POcorr,dIdPO = GetPrefOriDerv(refl,uniq,G,g,phfx,hfx,SGData,calcControls,parmDict)
    1099     for iPO in dIdPO:
    1100         dIdPO[iPO] /= POcorr
    11011174    dFdODF = {}
    11021175    dFdSA = [0,0,0]
    … …  
    11071180        for i in range(3):
    11081181            dFdSA[i] /= odfCor
    1109     dFdAb = GetAbsorbDerv(refl,hfx,calcControls,parmDict)
    1110     return dIdsh,dIdsp,dIdPola,dIdPO,dFdODF,dFdSA,dFdAb
     1182    elif calcControls[phfx+'poType'] == 'MD' or calcControls[phfx+'SHord']:
     1183        POcorr,dIdPO = GetPrefOriDerv(refl,uniq,G,g,phfx,hfx,SGData,calcControls,parmDict)       
     1184        for iPO in dIdPO:
     1185            dIdPO[iPO] /= POcorr
     1186    AbsCorr = GetAbsorb(refl,hfx,calcControls,parmDict)
     1187    dFdAb = GetAbsorbDerv(refl,hfx,calcControls,parmDict)/AbsCorr
     1188    dFdEx = GetPwdrExtDerv(refl,pfx,phfx,hfx,calcControls,parmDict)
     1189    return dIdsh,dIdsp,dIdPola,dIdPO,dFdODF,dFdSA,dFdAb,dFdEx
    11111190       
    11121191def GetSampleSigGam(refl,wave,G,GB,hfx,phfx,calcControls,parmDict):
    … …  
    15881667                refl[5] += GetHStrainShift(refl,SGData,phfx,parmDict)               #apply hydrostatic strain shift
    15891668                refl[6:8] = GetReflSigGamCW(refl,wave,G,GB,phfx,calcControls,parmDict)    #peak sig & gam
    1590                 refl[11] = GetIntensityCorr(refl,Uniq,G,g,pfx,phfx,hfx,SGData,calcControls,parmDict)*Vst*Lorenz
     1669                refl[11:15] = GetIntensityCorr(refl,Uniq,G,g,pfx,phfx,hfx,SGData,calcControls,parmDict)
     1670                refl[11] *= Vst*Lorenz
    15911671                 
    15921672                if Phase['General'].get('doPawley'):
    … …  
    16231703                refl[6:8] = GetReflSigGamTOF(refl,G,GB,phfx,calcControls,parmDict)    #peak sig & gam
    16241704                refl[12:14] = GetReflAlpBet(refl,hfx,parmDict)
    1625                 refl[11] = GetIntensityCorr(refl,Uniq,G,g,pfx,phfx,hfx,SGData,calcControls,parmDict)*Vst*Lorenz
     1705                refl[11],refl[15],refl[16],refl[17] = GetIntensityCorr(refl,Uniq,G,g,pfx,phfx,hfx,SGData,calcControls,parmDict)
     1706                refl[11] *= Vst*Lorenz
    16261707                if Phase['General'].get('doPawley'):
    16271708                    try:
    … …  
    17291810            h,k,l = refl[:3]
    17301811            Uniq = np.inner(refl[:3],SGMT)
    1731             dIdsh,dIdsp,dIdpola,dIdPO,dFdODF,dFdSA,dFdAb = GetIntensityDerv(refl,Uniq,G,g,pfx,phfx,hfx,SGData,calcControls,parmDict)
     1812            dIdsh,dIdsp,dIdpola,dIdPO,dFdODF,dFdSA,dFdAb,dFdEx = GetIntensityDerv(refl,Uniq,G,g,pfx,phfx,hfx,SGData,calcControls,parmDict)
    17321813            if 'C' in calcControls[hfx+'histType']:        #CW powder
    17331814                Wd,fmin,fmax = G2pwd.getWidthsCW(refl[5],refl[6],refl[7],shl)
    … …  
    17661847                dMdpk = np.zeros(shape=(6,lenBF))
    17671848                dMdipk = G2pwd.getdEpsVoigt(refl[5],refl[12],refl[13],refl[6],refl[7],ma.getdata(x[iBeg:iFin]))
    1768                 for i in range(5):
    1769                     dMdpk[i] += 100.*cw[iBeg:iFin]*refl[11]*refl[9]*dMdipk[i]
     1849                for i in range(6):
     1850                    dMdpk[i] += cw[iBeg:iFin]*refl[11]*refl[9]*dMdipk[i]
    17701851                dervDict = {'int':dMdpk[0],'pos':dMdpk[1],'alp':dMdpk[2],'bet':dMdpk[3],'sig':dMdpk[4],'gam':dMdpk[5]}           
    17711852            if Phase['General'].get('doPawley'):
    … …  
    18001881                    hfx+'alpha':[1./refl[4],'alp'],hfx+'beta-0':[1.0,'bet'],hfx+'beta-1':[1./refl[4]**4,'bet'],
    18011882                    hfx+'beta-q':[1./refl[4],'bet'],hfx+'sig-0':[1.0,'sig'],hfx+'sig-1':[refl[4]**2,'sig'],
    1802                     hfx+'sig-q':[refl[4],'sig'],hfx+'Absorption':[dFdAb,'int'],}
     1883                    hfx+'sig-q':[refl[4],'sig'],hfx+'Absorption':[dFdAb,'int'],hfx+'Extinction':[dFdEx,'int'],}
    18031884            for name in names:
    18041885                item = names[name]
    … …  
    19452026        for iref,ref in enumerate(refDict['RefList']):
    19462027            if ref[6] > 0:
    1947                 dervCor,dervDict = SCExtinction(ref,phfx,hfx,pfx,calcControls,parmDict,varylist) #puts correction in refl[11]
     2028                dervCor,dervDict = SCExtinction(ref,phfx,hfx,pfx,calcControls,parmDict,varylist)[1:]
    19482029                w = 1.0/ref[6]
    19492030                if w*ref[5] >= calcControls['minF/sig']:
    … …  
    19722053        for iref,ref in enumerate(refDict['RefList']):
    19732054            if ref[5] > 0.:
    1974                 dervCor,dervDict = SCExtinction(ref,phfx,hfx,pfx,calcControls,parmDict,varylist) #puts correction in refl[11]
     2055                dervCor,dervDict = SCExtinction(ref,phfx,hfx,pfx,calcControls,parmDict,varylist)[1:]
    19752056                Fo = np.sqrt(ref[5])
    19762057                Fc = np.sqrt(ref[7])
    … …  
    19802061                    for j,var in enumerate(varylist):
    19812062                        if var in dFdvDict:
    1982                             dMdvh[j][iref] = w*dFdvDict[var][iref]*parmDict[phfx+'Scale']*dervCor*ref[11]
     2063                            dMdvh[j][iref] = w*dFdvDict[var][iref]*parmDict[phfx+'Scale']*ref[11]       #*dervCor
    19832064                    for var in dependentVars:
    19842065                        if var in dFdvDict:
    1985                             depDerivDict[var][iref] = w*dFdvDict[var][iref]*parmDict[phfx+'Scale']*dervCor*ref[11]
     2066                            depDerivDict[var][iref] = w*dFdvDict[var][iref]*parmDict[phfx+'Scale']*ref[11]      #*dervCor
    19862067                    if phfx+'Scale' in varylist:
    1987                         dMdvh[varylist.index(phfx+'Scale')][iref] = w*ref[9]*dervCor
     2068                        dMdvh[varylist.index(phfx+'Scale')][iref] = w*ref[9]*ref[11]    #*dervCor
    19882069                    elif phfx+'Scale' in dependentVars:
    1989                         depDerivDict[phfx+'Scale'][iref] = w*ref[9]*dervCor                           
     2070                        depDerivDict[phfx+'Scale'][iref] = w*ref[9]*ref[11] #*dervCor                           
    19902071                    for item in ['Ep','Es','Eg']:
    19912072                        if phfx+item in varylist and dervDict:
    … …  
    22112292                for i,ref in enumerate(refDict['RefList']):
    22122293                    if ref[6] > 0:
    2213                         SCExtinction(ref,phfx,hfx,pfx,calcControls,parmDict,varylist) #puts correction in refl[11]
     2294                        ref[11] = SCExtinction(ref,phfx,hfx,pfx,calcControls,parmDict,varylist)[0]
    22142295                        w = 1.0/ref[6]
    22152296                        ref[7] = parmDict[phfx+'Scale']*ref[9]*ref[11]  #correct Fc^2 for extinction
    … …  
    22272308                for i,ref in enumerate(refDict['RefList']):
    22282309                    if ref[5] > 0.:
    2229                         SCExtinction(ref,phfx,hfx,pfx,calcControls,parmDict,varylist) #puts correction in refl[11]
    2230                         ref[7] = parmDict[phfx+'Scale']*ref[9]*ref[11]    #correct Fc^2 for extinctio
     2310                        ref[11] = SCExtinction(ref,phfx,hfx,pfx,calcControls,parmDict,varylist)[0]
     2311                        ref[7] = parmDict[phfx+'Scale']*ref[9]*ref[11]    #correct Fc^2 for extinction
    22312312                        ref[8] = ref[5]/(parmDict[phfx+'Scale']*ref[11])
    22322313                        Fo = np.sqrt(ref[5])
Note: See TracChangeset for help on using the changeset viewer.