Changeset 375


Ignore:
Timestamp:
Sep 19, 2011 1:20:32 PM (10 years ago)
Author:
vondreele
Message:

GSASIIgrid.py - add Pawley estimator to menu
GSASIIlattice.py - fix error in CosSinAngle?
GSASIIphsGUI.py - add Pawley estimator & modify size mustrain defaults
GSASIIpwd.py - modify polarization deriv
GSASIIstruct.py - modify refl index issue

Location:
trunk
Files:
5 edited

Legend:

Unmodified
Added
Removed
  • trunk/GSASIIgrid.py

    r373 r375  
    3737] = [wx.NewId() for _init_coll_MASK_Items in range(3)]
    3838
    39 [ wxID_PAWLEYLOAD, wxID_PAWLEYIMPORT, wxID_PAWLEYDELETE,
    40 ] = [wx.NewId() for _init_coll_PAWLEY_Items in range(3)]
     39[ wxID_PAWLEYLOAD, wxID_PAWLEYIMPORT, wxID_PAWLEYDELETE, wxID_PAWLEYESTIMATE,
     40] = [wx.NewId() for _init_coll_PAWLEY_Items in range(4)]
    4141
    4242[ wxID_INSTPRMRESET,wxID_CHANGEWAVETYPE,
     
    166166#        parent.Append(id=wxID_PAWLEYIMPORT, kind=wx.ITEM_NORMAL,text='Pawley import',
    167167#            help='Import Pawley reflection list')
     168        parent.Append(id=wxID_PAWLEYESTIMATE, kind=wx.ITEM_NORMAL,text='Pawley estimate',
     169            help='Estimate initial Pawley intensities')
    168170        parent.Append(id=wxID_PAWLEYDELETE, kind=wx.ITEM_NORMAL,text='Pawley delete',
    169171            help='Delete Pawley reflection list')
  • trunk/GSASIIlattice.py

    r360 r375  
    237237        cos(phi) & sin(phi)
    238238    """
    239     u = U/nl.norm(U)
    240     v = V/nl.norm(V)
     239    u = U/np.sqrt(np.inner(U,np.inner(G,U)))
     240    v = V/np.sqrt(np.inner(V,np.inner(G,V)))
    241241    cosP = np.inner(u,np.inner(G,v))
    242     sinP = np.sqrt(1.0-cosP**2)
     242    sinP = np.sqrt(max(0.0,1.0-cosP**2))
    243243    return cosP,sinP
    244244   
  • trunk/GSASIIphsGUI.py

    r372 r375  
    26272627                        UseList[histoName] = {'Histogram':histoName,'Show':False,
    26282628                            'Scale':[1.0,False],'Pref.Ori.':['MD',1.0,False,[0,0,1],0,{}],
    2629                             'Size':['isotropic',[10.,10.,],[False,False],[0,0,1],6*[0.0,],6*[False,]],
    2630                             'Mustrain':['isotropic',[1.0,1.0],[False,False],[0,0,1],
     2629                            'Size':['isotropic',[4.,4.,],[False,False],[0,0,1],6*[0.0,],6*[False,]],
     2630                            'Mustrain':['isotropic',[1000.0,1000.0],[False,False],[0,0,1],
    26312631                                NShkl*[0.01,],NShkl*[False,]],
    26322632                            'HStrain':[NDij*[0.0,],NDij*[False,]],                         
     
    27122712            wx.EndBusyCursor()
    27132713        data['Pawley ref'] = PawleyPeaks
     2714        FillPawleyReflectionsGrid()
     2715       
     2716    def OnPawleyEstimate(event):
     2717        try:
     2718            Refs = data['Pawley ref']
     2719            Histograms = data['Histograms']
     2720        except KeyError:
     2721            print '**** Error - no histograms defined for this phase ****'
     2722            return
     2723        HistoNames = Histograms.keys()
     2724        PatternId = G2gd.GetPatternTreeItemId(self,self.root,HistoNames[0])
     2725        xdata = self.PatternTree.GetItemPyData(PatternId)[1]
     2726        Inst = self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,PatternId,'Instrument Parameters'))
     2727        Inst = dict(zip(Inst[3],Inst[1]))
     2728        Sample = self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,PatternId,'Sample Parameters'))
     2729        if 'Lam' in Inst:
     2730            wave = Inst['Lam']
     2731        else:
     2732            wave = Inst['Lam1']
     2733       
     2734        posCorr = Inst['Zero']
     2735        print Sample
     2736        const = 9.e-2/(np.pi*Sample['Gonio. radius'])                  #shifts in microns
     2737       
     2738        for ref in Refs:
     2739            pos = 2.0*asind(wave/(2.0*ref[4]))
     2740            if 'Bragg' in Sample['Type']:
     2741                pos -= const*(4.*Sample['Shift'][0]*cosd(pos/2.0)+ \
     2742                    Sample['Transparency'][0]*sind(pos)*100.0)            #trans(=1/mueff) in cm
     2743            else:               #Debye-Scherrer - simple but maybe not right
     2744                pos -= const*(Sample['DisplaceX'][0]*cosd(pos)+Sample['DisplaceY'][0]*sind(pos))
     2745            indx = np.searchsorted(xdata[0],pos)
     2746            try:
     2747                ref[6] = xdata[1][indx]
     2748                print ref[:7],indx,pos
     2749            except IndexError:
     2750                pass
    27142751        FillPawleyReflectionsGrid()
    27152752                           
     
    27742811            self.dataFrame.SetMenuBar(self.dataFrame.PawleyMenu)
    27752812            self.dataFrame.Bind(wx.EVT_MENU, OnPawleyLoad, id=G2gd.wxID_PAWLEYLOAD)
     2813            self.dataFrame.Bind(wx.EVT_MENU, OnPawleyEstimate, id=G2gd.wxID_PAWLEYESTIMATE)
    27762814            self.dataFrame.Bind(wx.EVT_MENU, OnPawleyDelete, id=G2gd.wxID_PAWLEYDELETE)           
    27772815            FillPawleyReflectionsGrid()
  • trunk/GSASIIpwd.py

    r373 r375  
    137137    pola = ((1.0-Pola)*npcosd(Azm)**2+Pola*npsind(Azm)**2)*npcosd(Tth)**2+   \
    138138        (1.0-Pola)*npsind(Azm)**2+Pola*npcosd(Azm)**2
    139     dpdPola = npsind(Azm)**2*npsind(Tth)**2
    140     return pola,dpdPola
     139    dpdPola = -npsind(Tth)**2*(npsind(Azm)**2-npcosd(Azm)**2)
     140    return pola,dpdPola/pola
    141141   
    142142def Oblique(ObCoeff,Tth):
  • trunk/GSASIIstruct.py

    r374 r375  
    323323                if isum == jsum:
    324324                    G2mv.StoreEquivalence(varyJ,(varyI,))
     325                   
    325326def cellVary(pfx,SGData):
    326327    if SGData['SGLaue'] in ['-1',]:
     
    409410                    pawleyVary.append(pfx+'PWLref:'+str(i))
    410411            GetPawleyConstr(SGData['SGLaue'],PawleyRef,pawleyVary)      #does G2mv.StoreEquivalence
    411             phaseVary += pawleyVary   
     412            phaseVary += pawleyVary
    412413               
    413414    return phaseVary,phaseDict,pawleyLookup
     
    11581159        return gam
    11591160       
    1160     def GetIntensityCorr(refl,phfx,hfx,calcControls,parmDict):
     1161    def GetMarchDollase(refl,G,phfx,calcControls,parmDict):
     1162        MD = parmDict[phfx+'MD']
     1163        MDAxis = calcControls[phfx+'MDAxis']
     1164        sumMD = 0
     1165        for H in refl[10]:           
     1166            cosP,sinP = G2lat.CosSinAngle(H,MDAxis,G)
     1167            sumMD += 1/np.sqrt((MD*cosP)**2+sinP**2/MD)**3
     1168        return sumMD/len(refl[10])
     1169       
     1170    def GetIntensityCorr(refl,G,phfx,hfx,calcControls,parmDict):
    11611171        Icorr = parmDict[phfx+'Scale']*parmDict[hfx+'Scale']*refl[3]               #scale*multiplicity
    11621172        Icorr *= G2pwd.Polarization(parmDict[hfx+'Polariz.'],refl[5],parmDict[hfx+'Azimuth'])[0]
     1173        if calcControls[phfx+'poType'] == 'MD':
     1174            Icorr *= GetMarchDollase(refl,G,phfx,calcControls,parmDict)
    11631175       
    11641176        return Icorr
     
    12551267                refl[5] += GetHStrainShift(refl,SGData,phfx,parmDict)               #apply hydrostatic strain shift
    12561268                refl[6:8] = GetReflSIgGam(refl,wave,G,hfx,phfx,calcControls,parmDict,sizeEllipse)    #peak sig & gam
    1257                 Icorr = GetIntensityCorr(refl,phfx,hfx,calcControls,parmDict)
     1269                Icorr = GetIntensityCorr(refl,G,phfx,hfx,calcControls,parmDict)
    12581270                if 'Pawley' in Phase['General']['Type']:
    12591271                    try:
     
    12941306        #crystallite size derivatives
    12951307        if calcControls[phfx+'SizeType'] == 'isotropic':
    1296             gamDict[phfx+'Size:0'] = 1.80*wave/(np.pi*costh)
     1308            gamDict[phfx+'Size:0'] = -1.80*wave/(np.pi*costh)
    12971309        elif calcControls[phfx+'SizeType'] == 'uniaxial':
    12981310            H = np.array(refl[:3])
     
    13301342        return gamDict
    13311343       
    1332     def GetIntensityDerv(refl,phfx,hfx,calcControls,parmDict):
     1344    def GetMarchDollaseDerv(refl,G,phfx,calcControls,parmDict):
     1345        MD = parmDict[phfx+'MD']
     1346        MDAxis = calcControls[phfx+'MDAxis']
     1347        sumMD = 0
     1348        sumdMD = 0
     1349        for H in refl[10]:           
     1350            cosP,sinP = G2lat.CosSinAngle(H,MDAxis,G)
     1351            A = 1.0/np.sqrt((MD*cosP)**2+sinP**2/MD)
     1352            sumMD += A**3
     1353            sumdMD -= (1.5*A**5)*(2.0*MD*cosP**2-(sinP/MD)**2)
     1354        return sumMD/len(refl[10]),refl[8]*sumdMD
     1355       
     1356    def GetIntensityDerv(refl,G,phfx,hfx,calcControls,parmDict):
    13331357        Icorr = parmDict[phfx+'Scale']*parmDict[hfx+'Scale']*refl[3]               #scale*multiplicity
    13341358        pola,dpdPola = G2pwd.Polarization(parmDict[hfx+'Polariz.'],refl[5],parmDict[hfx+'Azimuth'])
    1335         dIdpola = Icorr*dpdPola
     1359        dIdpola = refl[8]*Icorr*dpdPola/refl[3]
    13361360        Icorr *= pola
     1361        MDcorr,dIdMD = GetMarchDollaseDerv(refl,G,phfx,calcControls,parmDict)
     1362        Icorr *= MDcorr
    13371363        dIdsh = Icorr/parmDict[hfx+'Scale']
    13381364        dIdsp = Icorr/parmDict[phfx+'Scale']
    1339        
    1340         return Icorr,dIdsh,dIdsp,dIdpola
     1365        dIdMD *= Icorr
     1366        return Icorr,dIdsh,dIdsp,dIdpola,dIdMD
    13411367       
    13421368    def GetReflPosDerv(refl,wave,A,hfx,calcControls,parmDict):
     
    14491475        if calcControls[phfx+'SizeType'] == 'ellipsoidal':
    14501476            sizeEllipse = G2lat.U6toUij([parmDIct[phfx+'Size:%d'%(i)] for i in range(6)])
    1451         for iref,refl in enumerate(refList):
     1477        for refl in refList:
    14521478            if 'C' in calcControls[hfx+'histType']:
    1453                 Icorr,dIdsh,dIdsp,dIdpola = GetIntensityDerv(refl,phfx,hfx,calcControls,parmDict)
     1479                h,k,l = refl[:3]
     1480                iref = pawleyLookup[pfx+'%d,%d,%d'%(h,k,l)]
     1481                Icorr,dIdsh,dIdsp,dIdpola,dIdMD = GetIntensityDerv(refl,G,phfx,hfx,calcControls,parmDict)
    14541482                hkl = refl[:3]
    14551483                pos = refl[5]
     
    14931521                    hfx+'I(L2)/I(L1)':[1.0,'L1/L2'],hfx+'Zero':[dpdZ,'pos'],hfx+'Lam':[dpdw,'pos'],
    14941522                    hfx+'Shift':[dpdSh,'pos'],hfx+'Transparency':[dpdTr,'pos'],hfx+'DisplaceX':[dpdX,'pos'],
    1495                     hfx+'DisplaceY':[dpdY,'pos'],}
     1523                    hfx+'DisplaceY':[dpdY,'pos'],phfx+'MD':[dIdMD,'int'],}
    14961524                for name in names:
    14971525                    if name in varylist:
     
    16441672#            for item in table: print item,table[item]               #useful debug - are things shifting?
    16451673            break                   #refinement succeeded - finish up!
    1646         except ValueError:          #result[1] is None on singular matrix
     1674        except TypeError:          #result[1] is None on singular matrix
    16471675            print '**** Refinement failed - singular matrix ****'
    16481676            Ipvt = result[2]['ipvt']
Note: See TracChangeset for help on using the changeset viewer.