Changeset 508


Ignore:
Timestamp:
Mar 2, 2012 3:40:56 PM (10 years ago)
Author:
vondreele
Message:

some refactoring in GSASIIphsGUI.py
implement LG mixing coeff for size & mustrain
remove MP stuff (it'll go in later)

Location:
trunk
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • trunk/GSASIIphsGUI.py

    r486 r508  
    253253        if 'Map' not in generalData:
    254254            generalData['Map'] = {'MapType':'','RefList':'','Resolution':4.0}
     255#        if 'SH Texture' not in generalData:
     256#            generalData['SH Texture'] = data['SH Texture']
    255257        generalData['NoAtoms'] = {}
    256258        generalData['BondRadii'] = []
     
    279281                    generalData['AtomMass'].append(Info['Isotopes'][generalData['Isotope'][atom[ct]]][0])
    280282                else:
     283                    generalData['Isotope'][atom[ct]] = 'Nat. Abund.'
    281284                    generalData['AtomMass'].append(Info['Mass'])
    282285                generalData['NoAtoms'][atom[ct]] = atom[cs-1]*float(atom[cs+1])
    283286                generalData['Color'].append(Info['Color'])
     287
     288################################################################################
     289##### General phase routines
     290################################################################################
    284291
    285292    def UpdateGeneral():
     
    516523            for elem in generalData['AtomTypes']:
    517524                choices = generalData['Isotopes'][elem].keys()
    518                 if elem not in generalData['Isotope']:
    519                     generalData['Isotope'][elem] = 'Nat. Abund.'
    520525                isoSel = wx.ComboBox(dataDisplay,-1,value=generalData['Isotope'][elem],choices=choices,
    521526                    style=wx.CB_READONLY|wx.CB_DROPDOWN)
     
    660665        dataDisplay.SetSize(Size)
    661666        G2frame.dataFrame.setSizePosLeft(Size)
     667
     668################################################################################
     669#####  Atom routines
     670################################################################################
    662671
    663672    def FillAtomsGrid():
     
    734743                            if parms == atomData[row][c]:
    735744                                Atoms.SelectRow(row,True)
     745                    SetupGeneral()
    736746                elif Atoms.GetColLabelValue(c) == 'residue':
    737747                    choice = []
     
    896906                if 'Atoms' in data['Drawing']:
    897907                    DrawAtomsReplaceByID(data['Drawing'],atomData[r],ID)
     908                SetupGeneral()
    898909            else:
    899910                event.Skip()
     
    10561067                FillAtomsGrid()
    10571068                G2plt.PlotStructure(G2frame,data)
     1069            SetupGeneral()
    10581070        event.StopPropagation()
    10591071
     
    12551267        for ind in indx:
    12561268            atomData[ind] = MakeDrawAtom(atom,atomData[ind])
     1269
     1270################################################################################
     1271##### Atom draw routines
     1272################################################################################
    12571273           
    12581274    def UpdateDrawAtoms():
     
    21432159        dataDisplay.SetSize(Size)
    21442160        G2frame.dataFrame.setSizePosLeft(Size)
     2161
     2162################################################################################
     2163####  Texture routines
     2164################################################################################
    21452165       
    21462166    def UpdateTexture():
     
    23692389        Size[1] = min(Size[1],450)
    23702390        G2frame.dataFrame.setSizePosLeft(Size)
     2391
     2392################################################################################
     2393##### DData routines
     2394################################################################################
    23712395       
    23722396    def UpdateDData():
     
    23842408        Indx = {}
    23852409       
    2386         def OnPlotSel(event):
    2387             Obj = event.GetEventObject()
    2388             generalData['Data plot type'] = Obj.GetStringSelection()
    2389             wx.CallAfter(UpdateDData)
    2390             G2plt.PlotSizeStrainPO(G2frame,data)
    2391            
    2392         def OnPOhkl(event):
    2393             Obj = event.GetEventObject()
    2394             Saxis = Obj.GetValue().split()
    2395             try:
    2396                 hkl = [int(Saxis[i]) for i in range(3)]
    2397             except (ValueError,IndexError):
    2398                 hkl = generalData['POhkl']
    2399             if not np.any(np.array(hkl)):
    2400                 hkl = generalData['POhkl']
    2401             generalData['POhkl'] = hkl
    2402             h,k,l = hkl
    2403             Obj.SetValue('%3d %3d %3d'%(h,k,l))
    2404             G2plt.PlotSizeStrainPO(G2frame,data)
    2405        
     2410        def PlotSizer():
     2411
     2412            def OnPlotSel(event):
     2413                Obj = event.GetEventObject()
     2414                generalData['Data plot type'] = Obj.GetStringSelection()
     2415                wx.CallAfter(UpdateDData)
     2416                G2plt.PlotSizeStrainPO(G2frame,data)
     2417               
     2418            def OnPOhkl(event):
     2419                Obj = event.GetEventObject()
     2420                Saxis = Obj.GetValue().split()
     2421                try:
     2422                    hkl = [int(Saxis[i]) for i in range(3)]
     2423                except (ValueError,IndexError):
     2424                    hkl = generalData['POhkl']
     2425                if not np.any(np.array(hkl)):
     2426                    hkl = generalData['POhkl']
     2427                generalData['POhkl'] = hkl
     2428                h,k,l = hkl
     2429                Obj.SetValue('%3d %3d %3d'%(h,k,l))
     2430                G2plt.PlotSizeStrainPO(G2frame,data)
     2431           
     2432            plotSizer = wx.BoxSizer(wx.VERTICAL)
     2433            choice = ['Mustrain','Size','Preferred orientation']
     2434            plotSel = wx.RadioBox(DData,-1,'Select plot type:',choices=choice,
     2435                majorDimension=3,style=wx.RA_SPECIFY_COLS)
     2436            plotSel.SetStringSelection(generalData['Data plot type'])
     2437            plotSel.Bind(wx.EVT_RADIOBOX,OnPlotSel)   
     2438            plotSizer.Add(plotSel)
     2439            if generalData['Data plot type'] == 'Preferred orientation':
     2440                POhklSizer = wx.BoxSizer(wx.HORIZONTAL)
     2441                POhklSizer.Add(wx.StaticText(DData,-1,' Plot preferred orientation for H K L: '),0,wx.ALIGN_CENTER_VERTICAL)
     2442                h,k,l = generalData['POhkl']
     2443                poAxis = wx.TextCtrl(DData,-1,'%3d %3d %3d'%(h,k,l),style=wx.TE_PROCESS_ENTER)
     2444                poAxis.Bind(wx.EVT_TEXT_ENTER,OnPOhkl)
     2445                poAxis.Bind(wx.EVT_KILL_FOCUS,OnPOhkl)
     2446                POhklSizer.Add(poAxis,0,wx.ALIGN_CENTER_VERTICAL)
     2447                plotSizer.Add(POhklSizer)           
     2448            return plotSizer
     2449           
     2450        def ScaleSizer():
     2451           
     2452            def OnScaleRef(event):
     2453                Obj = event.GetEventObject()
     2454                UseList[Indx[Obj.GetId()]]['Scale'][1] = Obj.GetValue()
     2455               
     2456            def OnScaleVal(event):
     2457                Obj = event.GetEventObject()
     2458                try:
     2459                    scale = float(Obj.GetValue())
     2460                    if scale > 0:
     2461                        UseList[Indx[Obj.GetId()]]['Scale'][0] = scale
     2462                except ValueError:
     2463                    pass
     2464                Obj.SetValue("%.4f"%(UseList[Indx[Obj.GetId()]]['Scale'][0]))          #reset in case of error
     2465                           
     2466            scaleSizer = wx.BoxSizer(wx.HORIZONTAL)
     2467            scaleRef = wx.CheckBox(DData,-1,label=' Phase fraction: ')
     2468            scaleRef.SetValue(UseList[item]['Scale'][1])
     2469            Indx[scaleRef.GetId()] = item
     2470            scaleRef.Bind(wx.EVT_CHECKBOX, OnScaleRef)
     2471            scaleSizer.Add(scaleRef,0,wx.ALIGN_CENTER_VERTICAL)
     2472            scaleVal = wx.TextCtrl(DData,wx.ID_ANY,
     2473                '%.4f'%(UseList[item]['Scale'][0]),style=wx.TE_PROCESS_ENTER)
     2474            Indx[scaleVal.GetId()] = item
     2475            scaleVal.Bind(wx.EVT_TEXT_ENTER,OnScaleVal)
     2476            scaleVal.Bind(wx.EVT_KILL_FOCUS,OnScaleVal)
     2477            scaleSizer.Add(scaleVal,0,wx.ALIGN_CENTER_VERTICAL)
     2478            return scaleSizer
     2479           
    24062480        def OnShowData(event):
    24072481            Obj = event.GetEventObject()
     
    24392513                    dlg.Destroy()
    24402514           
    2441         def OnScaleRef(event):
     2515        def OnLGmixRef(event):
    24422516            Obj = event.GetEventObject()
    2443             UseList[Indx[Obj.GetId()]]['Scale'][1] = Obj.GetValue()
    2444            
    2445         def OnScaleVal(event):
     2517            hist,name = Indx[Obj.GetId()]
     2518            UseList[hist][name][2][2] = Obj.GetValue()
     2519           
     2520        def OnLGmixVal(event):
    24462521            Obj = event.GetEventObject()
     2522            hist,name = Indx[Obj.GetId()]
    24472523            try:
    2448                 scale = float(Obj.GetValue())
    2449                 if scale > 0:
    2450                     UseList[Indx[Obj.GetId()]]['Scale'][0] = scale
     2524                value = float(Obj.GetValue())
     2525                if 0.1 <= value <= 1:
     2526                    UseList[hist][name][1][2] = value
     2527                else:
     2528                    raise ValueError
    24512529            except ValueError:
    24522530                pass
    2453             Obj.SetValue("%.4f"%(UseList[Indx[Obj.GetId()]]['Scale'][0]))          #reset in case of error
    2454            
     2531            Obj.SetValue("%.4f"%(UseList[hist][name][1][2]))          #reset in case of error
     2532
    24552533        def OnSizeType(event):
    24562534            Obj = event.GetEventObject()
     
    25092587            hist = Indx[Obj.GetId()]
    25102588            UseList[hist]['Mustrain'][0] = Obj.GetValue()
     2589            wx.CallAfter(UpdateDData)
    25112590            G2plt.PlotSizeStrainPO(G2frame,data)
    2512             wx.CallAfter(UpdateDData)
    25132591           
    25142592        def OnStrainRef(event):
     
    26482726            return axis
    26492727           
    2650         def PlotSizer():
    2651             plotSizer = wx.BoxSizer(wx.VERTICAL)
    2652             choice = ['Mustrain','Size','Preferred orientation']
    2653             plotSel = wx.RadioBox(DData,-1,'Select plot type:',choices=choice,
    2654                 majorDimension=3,style=wx.RA_SPECIFY_COLS)
    2655             plotSel.SetStringSelection(generalData['Data plot type'])
    2656             plotSel.Bind(wx.EVT_RADIOBOX,OnPlotSel)   
    2657             plotSizer.Add(plotSel)
    2658             if generalData['Data plot type'] == 'Preferred orientation':
    2659                 POhklSizer = wx.BoxSizer(wx.HORIZONTAL)
    2660                 POhklSizer.Add(wx.StaticText(DData,-1,' Plot preferred orientation for H K L: '),0,wx.ALIGN_CENTER_VERTICAL)
    2661                 h,k,l = generalData['POhkl']
    2662                 poAxis = wx.TextCtrl(DData,-1,'%3d %3d %3d'%(h,k,l),style=wx.TE_PROCESS_ENTER)
    2663                 poAxis.Bind(wx.EVT_TEXT_ENTER,OnPOhkl)
    2664                 poAxis.Bind(wx.EVT_KILL_FOCUS,OnPOhkl)
    2665                 POhklSizer.Add(poAxis,0,wx.ALIGN_CENTER_VERTICAL)
    2666                 plotSizer.Add(POhklSizer)           
    2667             return plotSizer
    2668            
    2669         def ScaleSizer():
    2670             scaleSizer = wx.BoxSizer(wx.HORIZONTAL)
    2671             scaleRef = wx.CheckBox(DData,-1,label=' Phase fraction: ')
    2672             scaleRef.SetValue(UseList[item]['Scale'][1])
    2673             Indx[scaleRef.GetId()] = item
    2674             scaleRef.Bind(wx.EVT_CHECKBOX, OnScaleRef)
    2675             scaleSizer.Add(scaleRef,0,wx.ALIGN_CENTER_VERTICAL)
    2676             scaleVal = wx.TextCtrl(DData,wx.ID_ANY,
    2677                 '%.4f'%(UseList[item]['Scale'][0]),style=wx.TE_PROCESS_ENTER)
    2678             Indx[scaleVal.GetId()] = item
    2679             scaleVal.Bind(wx.EVT_TEXT_ENTER,OnScaleVal)
    2680             scaleVal.Bind(wx.EVT_KILL_FOCUS,OnScaleVal)
    2681             scaleSizer.Add(scaleVal,0,wx.ALIGN_CENTER_VERTICAL)
    2682             return scaleSizer
    2683            
    26842728        def TopSizer(name,choices,parm,OnType):
    26852729            topSizer = wx.BoxSizer(wx.HORIZONTAL)
     
    26932737            return topSizer
    26942738           
     2739        def LGmixSizer(name,OnVal,OnRef):
     2740            lgmixSizer = wx.BoxSizer(wx.HORIZONTAL)
     2741            lgmixRef = wx.CheckBox(DData,-1,label='LGmix')
     2742            lgmixRef.thisown = False
     2743            lgmixRef.SetValue(UseList[item][name][2][2])
     2744            Indx[lgmixRef.GetId()] = [item,name]
     2745            lgmixRef.Bind(wx.EVT_CHECKBOX, OnRef)
     2746            lgmixSizer.Add(lgmixRef,0,wx.ALIGN_CENTER_VERTICAL)
     2747            lgmixVal = wx.TextCtrl(DData,wx.ID_ANY,
     2748                '%.4f'%(UseList[item][name][1][2]),style=wx.TE_PROCESS_ENTER)
     2749            Indx[lgmixVal.GetId()] = [item,name]
     2750            lgmixVal.Bind(wx.EVT_TEXT_ENTER,OnVal)
     2751            lgmixVal.Bind(wx.EVT_KILL_FOCUS,OnVal)
     2752            lgmixSizer.Add(lgmixVal,0,wx.ALIGN_CENTER_VERTICAL)
     2753            return lgmixSizer
     2754                       
    26952755        def IsoSizer(name,parm,fmt,OnVal,OnRef):
    26962756            isoSizer = wx.BoxSizer(wx.HORIZONTAL)
     
    28832943        for item in keyList:
    28842944            histData = UseList[item]
    2885            
     2945###### Patch to add LGmix to Size & Mustrain
     2946            if len(histData['Size'][1]) == 2:
     2947                histData['Size'][1].append(0.6667)
     2948                histData['Size'][2].append(False)
     2949                histData['Mustrain'][1].append(0.6667)
     2950                histData['Mustrain'][2].append(False)
     2951                UseList[item] = histData
     2952###### end patch
    28862953            showSizer = wx.BoxSizer(wx.HORIZONTAL)
    28872954            showData = wx.CheckBox(DData,-1,label=' Show '+item)
     
    29092976                    isoSizer.Add(IsoSizer(u' Cryst. size(\xb5m): ','Size','%.3f',
    29102977                        OnSizeVal,OnSizeRef),0,wx.ALIGN_CENTER_VERTICAL)
     2978                    isoSizer.Add(LGmixSizer('Size',OnLGmixVal,OnLGmixRef))
    29112979                    mainSizer.Add(isoSizer)
    29122980                elif UseList[item]['Size'][0] == 'uniaxial':
     
    29152983                        'Size',OnSizeType),0,wx.ALIGN_CENTER_VERTICAL)
    29162984                    uniSizer.Add(UniSizer('Size',OnSizeAxis),0,wx.ALIGN_CENTER_VERTICAL)
     2985                    uniSizer.Add(LGmixSizer('Size',OnLGmixVal,OnLGmixRef))
    29172986                    mainSizer.Add(uniSizer)
    29182987                    mainSizer.Add(UniDataSizer(u'size(\xb5m): ','Size','%.3f',OnSizeVal,OnSizeRef))
     
    29222991                        'Size',OnSizeType),0,wx.ALIGN_CENTER_VERTICAL)
    29232992                    ellSizer.Add(wx.StaticText(DData,-1,u' Coefficients(\xb5m): '),0,wx.ALIGN_CENTER_VERTICAL)
     2993                    ellSizer.Add(LGmixSizer('Size',OnLGmixVal,OnLGmixRef))
    29242994                    mainSizer.Add(ellSizer)
    29252995                    mainSizer.Add(EllSizeDataSizer())
     
    29323002                    isoSizer.Add(IsoSizer(' microstrain: ','Mustrain','%.1f',
    29333003                        OnStrainVal,OnStrainRef),0,wx.ALIGN_CENTER_VERTICAL)                   
     3004                    isoSizer.Add(LGmixSizer('Mustrain',OnLGmixVal,OnLGmixRef))
    29343005                    mainSizer.Add(isoSizer)
    29353006                    mainSizer.Add((0,5),0)
     
    29393010                        'Mustrain',OnStrainType),0,wx.ALIGN_CENTER_VERTICAL)
    29403011                    uniSizer.Add(UniSizer('Mustrain',OnStrainAxis),0,wx.ALIGN_CENTER_VERTICAL)
     3012                    uniSizer.Add(LGmixSizer('Mustrain',OnLGmixVal,OnLGmixRef))
    29413013                    mainSizer.Add(uniSizer)
    29423014                    mainSizer.Add(UniDataSizer('mustrain: ','Mustrain','%.1f',OnStrainVal,OnStrainRef))
     
    29453017                    genSizer.Add(TopSizer(' Mustrain model: ',['isotropic','uniaxial','generalized',],
    29463018                        'Mustrain',OnStrainType),0,wx.ALIGN_CENTER_VERTICAL)
    2947                     genSizer.Add(wx.StaticText(DData,-1,' Coefficients: '),0,wx.ALIGN_CENTER_VERTICAL)
     3019                    genSizer.Add(LGmixSizer('Mustrain',OnLGmixVal,OnLGmixRef))
    29483020                    mainSizer.Add(genSizer)
    29493021                    mainSizer.Add(GenStrainDataSizer())                       
     
    30333105                        UseList[histoName] = {'Histogram':histoName,'Show':False,
    30343106                            'Scale':[1.0,False],'Pref.Ori.':['MD',1.0,False,[0,0,1],0,{}],
    3035                             'Size':['isotropic',[4.,4.,],[False,False],[0,0,1],[4.,4.,4.,0.,0.,0.],6*[False,]],
    3036                             'Mustrain':['isotropic',[1000.0,1000.0],[False,False],[0,0,1],
    3037                                 NShkl*[0.01,],NShkl*[False,]],
     3107                            'Size':['isotropic',[4.,4.,0.66667],[False,False,False],[0,0,1],
     3108                                [4.,4.,4.,0.,0.,0.],6*[False,]],
     3109                            'Mustrain':['isotropic',[1000.0,1000.0,0.666667],[False,False,False],[0,0,1],
     3110                                (NShkl)*[0.01,],NShkl*[False,]],
    30383111                            'HStrain':[NDij*[0.0,],NDij*[False,]],                         
    30393112                            'Extinction':[0.0,False]}
     
    30673140            finally:
    30683141                dlg.Destroy()
     3142
     3143################################################################################
     3144##### Pawley routines
     3145################################################################################
    30693146
    30703147    def FillPawleyReflectionsGrid():
     
    31743251            data['Pawley ref'] = []
    31753252            FillPawleyReflectionsGrid()
     3253
     3254################################################################################
     3255##### End of main routines
     3256################################################################################
    31763257   
    31773258    def OnFourierMaps(event):
     
    32723353        DData = wx.ScrolledWindow(G2frame.dataDisplay)
    32733354        G2frame.dataDisplay.AddPage(DData,'Data')
    3274         Texture = wx.ScrolledWindow(G2frame.dataDisplay)
    3275         G2frame.dataDisplay.AddPage(Texture,'Texture')
    32763355        Atoms = G2gd.GSGrid(G2frame.dataDisplay)
    32773356        G2frame.dataDisplay.AddPage(Atoms,'Atoms')
     
    32803359        drawAtoms = G2gd.GSGrid(G2frame.dataDisplay)
    32813360        G2frame.dataDisplay.AddPage(drawAtoms,'Draw Atoms')
     3361        Texture = wx.ScrolledWindow(G2frame.dataDisplay)
     3362        G2frame.dataDisplay.AddPage(Texture,'Texture')
    32823363
    32833364    G2frame.dataDisplay.Bind(wx.EVT_NOTEBOOK_PAGE_CHANGED, OnPageChanged)
  • trunk/GSASIIstruct.py

    r495 r508  
    1313import math
    1414import cPickle
    15 import multiprocessing as mp
    1615import numpy as np
    1716import numpy.linalg as nl
     
    3130acosd = lambda x: 180.*np.arccos(x)/np.pi
    3231atan2d = lambda y,x: 180.*np.arctan2(y,x)/np.pi
     32   
     33ateln2 = 8.0*math.log(2.0)
    3334
    3435
     
    404405    if 'Hessian' in Controls['deriv type']:
    405406        print ' Maximum number of cycles:',Controls['max cyc']
    406         print ' Maximum histogram processes:',Controls['max Hprocess']
    407         print ' Maximum reflection processes:',Controls['max Rprocess']
    408407    else:
    409408        print ' Minimum delta-M/M for convergence: ','%.2g'%(Controls['min dM/M'])
     
    642641#            elif General['Type'] == 'macromolecular':
    643642
    644                
    645643            textureData = General['SH Texture']
    646644            if textureData['Order']:
     
    965963            if hapData[0] == 'uniaxial':
    966964                line += ' axial:'+'%12.3f'%(hapData[1][1])+' Refine? '+str(hapData[2][1])
     965            line += '\n\t LG mixing coeff.: %12.4f'%(hapData[1][2])+' Refine? '+str(hapData[2][2])
    967966            print line
    968967        else:
    969             print '\n Size model    : %s'%(hapData[0])
     968            print '\n Size model    : %s'%(hapData[0])+ \
     969                '\n\t LG mixing coeff.:%12.4f'%(hapData[1][2])+' Refine? '+str(hapData[2][2])
    970970            Snames = ['S11','S22','S33','S12','S13','S23']
    971971            ptlbls = ' names :'
     
    986986            if hapData[0] == 'uniaxial':
    987987                line += ' axial:'+'%12.1f'%(hapData[1][1])+' Refine? '+str(hapData[2][1])
     988            line +='\n\t LG mixing coeff.:%12.4f'%(hapData[1][2])+' Refine? '+str(hapData[2][2])
    988989            print line
    989990        else:
    990             print '\n Mustrain model: %s'%(hapData[0])
     991            print '\n Mustrain model: %s'%(hapData[0])+ \
     992                '\n\t LG mixing coeff.:%12.4f'%(hapData[1][2])+' Refine? '+str(hapData[2][2])
    991993            Snames = G2spc.MustrainNames(SGData)
    992994            ptlbls = ' names :'
     
    10851087            for item in ['Mustrain','Size']:
    10861088                controlDict[pfx+item+'Type'] = hapData[item][0]
     1089                hapDict[pfx+item+':mx'] = hapData[item][1][2]
     1090                if hapData[item][2][2]:
     1091                    hapVary.append(pfx+item+':mx')
    10871092                if hapData[item][0] in ['isotropic','uniaxial']:
    10881093                    hapDict[pfx+item+':i'] = hapData[item][1][0]
     
    10951100                            hapVary.append(pfx+item+':a')
    10961101                else:       #generalized for mustrain or ellipsoidal for size
     1102                    Nterms = len(hapData[item][4])
    10971103                    if item == 'Mustrain':
    10981104                        names = G2spc.MustrainNames(SGData)
     
    11021108                            pwrs.append([int(h),int(k),int(l)])
    11031109                        controlDict[pfx+'MuPwrs'] = pwrs
    1104                     for i in range(len(hapData[item][4])):
     1110                    for i in range(Nterms):
    11051111                        sfx = ':'+str(i)
    11061112                        hapDict[pfx+item+sfx] = hapData[item][4][i]
     
    11461152            line += ' equatorial:%12.3f'%(hapData[1][0])
    11471153            if sizeSig[0][0]:
    1148                 line += ', sig: %8.3f'%(sizeSig[0][0])
     1154                line += ', sig:%8.3f'%(sizeSig[0][0])
    11491155                refine = True
    11501156            if hapData[0] == 'uniaxial':
     
    11521158                if sizeSig[0][1]:
    11531159                    refine = True
    1154                     line += ', sig: %8.3f'%(sizeSig[0][1])
    1155             if refine:
    1156                 print line
     1160                    line += ', sig:%8.3f'%(sizeSig[0][1])
     1161            line += ' LG mix coeff.:%12.4f'%(hapData[1][2])
     1162            if sizeSig[0][2]:
     1163                refine = True
     1164                line += ', sig:%8.3f'%(sizeSig[0][2])
    11571165        else:
     1166            line += ' LG mix coeff.:%12.4f'%(hapData[1][2])
     1167            if sizeSig[0][2]:
     1168                refine = True
     1169                line += ', sig:%8.3f'%(sizeSig[0][2])
    11581170            Snames = ['S11','S22','S33','S12','S13','S23']
    11591171            ptlbls = ' name  :'
     
    11851197                line += ' axial:%12.1f'%(hapData[1][1])
    11861198                if mustrainSig[0][1]:
    1187                      line += ', sig: %8.1f'%(mustrainSig[0][1])
    1188             if refine:
    1189                 print line
     1199                     line += ', sig:%8.1f'%(mustrainSig[0][1])
     1200            line += ' LG mix coeff.:%12.4f'%(hapData[1][2])
     1201            if mustrainSig[0][2]:
     1202                refine = True
     1203                line += ', sig:%8.3f'%(mustrainSig[0][2])
    11901204        else:
     1205            line += ' LG mix coeff.:%12.4f'%(hapData[1][2])
     1206            if mustrainSig[0][2]:
     1207                refine = True
     1208                line += ', sig:%8.3f'%(mustrainSig[0][2])
    11911209            Snames = G2spc.MustrainNames(SGData)
    11921210            ptlbls = ' name  :'
     
    12871305                else:
    12881306                    PrintSHPOAndSig(hapData['Pref.Ori.'],PhFrExtPOSig)
    1289             SizeMuStrSig = {'Mustrain':[[0,0],[0 for i in range(len(hapData['Mustrain'][4]))]],
    1290                 'Size':[[0,0],[0 for i in range(len(hapData['Size'][4]))]],
     1307            SizeMuStrSig = {'Mustrain':[[0,0,0],[0 for i in range(len(hapData['Mustrain'][4]))]],
     1308                'Size':[[0,0,0],[0 for i in range(len(hapData['Size'][4]))]],
    12911309                'HStrain':{}}                 
    12921310            for item in ['Mustrain','Size']:
     1311                hapData[item][1][2] = parmDict[pfx+item+':mx']
     1312                hapData[item][1][2] = min(1.,max(0.1,hapData[item][1][2]))
     1313                if pfx+item+':mx' in sigDict:
     1314                    SizeMuStrSig[item][0][2] = sigDict[pfx+item+':mx']
    12931315                if hapData[item][0] in ['isotropic','uniaxial']:                   
    12941316                    hapData[item][1][0] = parmDict[pfx+item+':i']
     
    13041326                            SizeMuStrSig[item][0][1] = sigDict[pfx+item+':a']
    13051327                else:       #generalized for mustrain or ellipsoidal for size
    1306                     for i in range(len(hapData[item][4])):
     1328                    Nterms = len(hapData[item][4])
     1329                    for i in range(Nterms):
    13071330                        sfx = ':'+str(i)
    13081331                        hapData[item][4][i] = parmDict[pfx+item+sfx]
     
    16811704    BLvals = {}
    16821705    for El in BLtables:
    1683         BLvals[El] = BLtables[El][1]
     1706        BLvals[El] = BLtables[El][1][1]
    16841707    return BLvals
    16851708       
     
    17201743            else:       #'X'
    17211744                refl[-1] = getFFvalues(FFtables,SQ)
    1722         for i,El in enumerate(Tdata):           
     1745        for i,El in enumerate(Tdata):
    17231746            FF[i] = refl[-1][El]           
    17241747        SQfactor = 4.0*SQ*twopisq
     
    20212044    return dIdsh,dIdsp,dIdPola,dIdPO,dFdODF,dFdSA
    20222045       
    2023 def GetSampleGam(refl,wave,G,GB,phfx,calcControls,parmDict):
     2046def GetSampleSigGam(refl,wave,G,GB,phfx,calcControls,parmDict):
    20242047    costh = cosd(refl[5]/2.)
    20252048    #crystallite size
    20262049    if calcControls[phfx+'SizeType'] == 'isotropic':
    2027         gam = 1.8*wave/(np.pi*parmDict[phfx+'Size:i']*costh)
     2050        Sgam = 1.8*wave/(np.pi*parmDict[phfx+'Size:i']*costh)
    20282051    elif calcControls[phfx+'SizeType'] == 'uniaxial':
    20292052        H = np.array(refl[:3])
    20302053        P = np.array(calcControls[phfx+'SizeAxis'])
    20312054        cosP,sinP = G2lat.CosSinAngle(H,P,G)
    2032         gam = (1.8*wave/np.pi)/(parmDict[phfx+'Size:i']*parmDict[phfx+'Size:a']*costh)
    2033         gam *= np.sqrt((sinP*parmDict[phfx+'Size:a'])**2+(cosP*parmDict[phfx+'Size:i'])**2)
     2055        Sgam = (1.8*wave/np.pi)/(parmDict[phfx+'Size:i']*parmDict[phfx+'Size:a']*costh)
     2056        Sgam *= np.sqrt((sinP*parmDict[phfx+'Size:a'])**2+(cosP*parmDict[phfx+'Size:i'])**2)
    20342057    else:           #ellipsoidal crystallites
    20352058        Sij =[parmDict[phfx+'Size:%d'%(i)] for i in range(6)]
    20362059        H = np.array(refl[:3])
    20372060        lenR = G2pwd.ellipseSize(H,Sij,GB)
    2038         gam = 1.8*wave/(np.pi*costh*lenR)
     2061        Sgam = 1.8*wave/(np.pi*costh*lenR)
    20392062    #microstrain               
    20402063    if calcControls[phfx+'MustrainType'] == 'isotropic':
    2041         gam += 0.018*parmDict[phfx+'Mustrain:i']*tand(refl[5]/2.)/np.pi
     2064        Mgam = 0.018*parmDict[phfx+'Mustrain:i']*tand(refl[5]/2.)/np.pi
    20422065    elif calcControls[phfx+'MustrainType'] == 'uniaxial':
    20432066        H = np.array(refl[:3])
     
    20462069        Si = parmDict[phfx+'Mustrain:i']
    20472070        Sa = parmDict[phfx+'Mustrain:a']
    2048         gam += 0.018*Si*Sa*tand(refl[5]/2.)/(np.pi*np.sqrt((Si*cosP)**2+(Sa*sinP)**2))
     2071        Mgam = 0.018*Si*Sa*tand(refl[5]/2.)/(np.pi*np.sqrt((Si*cosP)**2+(Sa*sinP)**2))
    20492072    else:       #generalized - P.W. Stephens model
    20502073        pwrs = calcControls[phfx+'MuPwrs']
     
    20522075        for i,pwr in enumerate(pwrs):
    20532076            sum += parmDict[phfx+'Mustrain:'+str(i)]*refl[0]**pwr[0]*refl[1]**pwr[1]*refl[2]**pwr[2]
    2054         gam += 0.018*refl[4]**2*tand(refl[5]/2.)*sum           
    2055     return gam
    2056        
    2057 def GetSampleGamDerv(refl,wave,G,GB,phfx,calcControls,parmDict):
     2077        Mgam = 0.018*refl[4]**2*tand(refl[5]/2.)*sum
     2078    gam = Sgam*parmDict[phfx+'Size:mx']+Mgam*parmDict[phfx+'Mustrain:mx']
     2079    sig = (Sgam*(1.-parmDict[phfx+'Size:mx']))**2+(Mgam*(1.-parmDict[phfx+'Mustrain:mx']))**2
     2080    sig /= ateln2
     2081    return sig,gam
     2082       
     2083def GetSampleSigGamDerv(refl,wave,G,GB,phfx,calcControls,parmDict):
    20582084    gamDict = {}
     2085    sigDict = {}
    20592086    costh = cosd(refl[5]/2.)
    20602087    tanth = tand(refl[5]/2.)
    20612088    #crystallite size derivatives
    20622089    if calcControls[phfx+'SizeType'] == 'isotropic':
    2063         gamDict[phfx+'Size:i'] = -1.80*wave/(np.pi*costh)
     2090        Sgam = 1.8*wave/(np.pi*parmDict[phfx+'Size:i']*costh)
     2091        gamDict[phfx+'Size:i'] = -1.80*wave*parmDict[phfx+'Size:mx']/(np.pi*costh)
     2092        sigDict[phfx+'Size:i'] = -3.60*Sgam*wave*(1.-parmDict[phfx+'Size:mx'])**2/(np.pi*costh*ateln2)
    20642093    elif calcControls[phfx+'SizeType'] == 'uniaxial':
    20652094        H = np.array(refl[:3])
     
    20702099        gami = (1.8*wave/np.pi)/(Si*Sa)
    20712100        sqtrm = np.sqrt((sinP*Sa)**2+(cosP*Si)**2)
    2072         gam = gami*sqtrm/costh           
    2073         gamDict[phfx+'Size:i'] = gami*Si*cosP**2/(sqtrm*costh)-gam/Si
    2074         gamDict[phfx+'Size:a'] = gami*Sa*sinP**2/(sqtrm*costh)-gam/Sa         
     2101        Sgam = gami*sqtrm
     2102        gam = Sgam/costh
     2103        dsi = (gami*Si*cosP**2/(sqtrm*costh)-gam/Si)
     2104        dsa = (gami*Sa*sinP**2/(sqtrm*costh)-gam/Sa)
     2105        gamDict[phfx+'Size:i'] = dsi*parmDict[phfx+'Size:mx']
     2106        gamDict[phfx+'Size:a'] = dsa*parmDict[phfx+'Size:mx']
     2107        sigDict[phfx+'Size:i'] = 2.*dsi*Sgam*(1.-parmDict[phfx+'Size:mx'])**2/ateln2
     2108        sigDict[phfx+'Size:a'] = 2.*dsa*Sgam*(1.-parmDict[phfx+'Size:mx'])**2/ateln2
    20752109    else:           #ellipsoidal crystallites
    20762110        const = 1.8*wave/(np.pi*costh)
    20772111        Sij =[parmDict[phfx+'Size:%d'%(i)] for i in range(6)]
    20782112        H = np.array(refl[:3])
    2079         R,dRdS = G2pwd.ellipseSizeDerv(H,Sij,GB)
     2113        lenR,dRdS = G2pwd.ellipseSizeDerv(H,Sij,GB)
     2114        Sgam = 1.8*wave/(np.pi*costh*lenR)
    20802115        for i,item in enumerate([phfx+'Size:%d'%(j) for j in range(6)]):
    2081             gamDict[item] = -(const/R**2)*dRdS[i]
     2116            gamDict[item] = -(const/lenR**2)*dRdS[i]*parmDict[phfx+'Size:mx']
     2117            sigDict[item] = -2.*Sgam*(const/lenR**2)*dRdS[i]*(1.-parmDict[phfx+'Size:mx'])**2/ateln2
     2118    gamDict[phfx+'Size:mx'] = Sgam
     2119    sigDict[phfx+'Size:mx'] = -2.*Sgam**2*(1.-parmDict[phfx+'Size:mx'])/ateln2
     2120           
    20822121    #microstrain derivatives               
    20832122    if calcControls[phfx+'MustrainType'] == 'isotropic':
    2084         gamDict[phfx+'Mustrain:i'] =  0.018*tanth/np.pi           
     2123        Mgam = 0.018*parmDict[phfx+'Mustrain:i']*tand(refl[5]/2.)/np.pi
     2124        gamDict[phfx+'Mustrain:i'] =  0.018*tanth*parmDict[phfx+'Mustrain:mx']/np.pi
     2125        sigDict[phfx+'Mustrain:i'] =  0.036*Mgam*tanth*(1.-parmDict[phfx+'Mustrain:mx'])**2/(np.pi*ateln2)       
    20852126    elif calcControls[phfx+'MustrainType'] == 'uniaxial':
    20862127        H = np.array(refl[:3])
     
    20912132        gami = 0.018*Si*Sa*tanth/np.pi
    20922133        sqtrm = np.sqrt((Si*cosP)**2+(Sa*sinP)**2)
    2093         gam = gami/sqtrm
    2094         gamDict[phfx+'Mustrain:i'] = gam/Si-gami*Si*cosP**2/sqtrm**3
    2095         gamDict[phfx+'Mustrain:a'] = gam/Sa-gami*Sa*sinP**2/sqtrm**3
     2134        Mgam = gami/sqtrm
     2135        dsi = -gami*Si*cosP**2/sqtrm**3
     2136        dsa = -gami*Sa*sinP**2/sqtrm**3
     2137        gamDict[phfx+'Mustrain:i'] = (Mgam/Si+dsi)*parmDict[phfx+'Mustrain:mx']
     2138        gamDict[phfx+'Mustrain:a'] = (Mgam/Sa+dsa)*parmDict[phfx+'Mustrain:mx']
     2139        sigDict[phfx+'Mustrain:i'] = 2*(Mgam/Si+dsi)*Mgam*(1.-parmDict[phfx+'Mustrain:mx'])**2/ateln2
     2140        sigDict[phfx+'Mustrain:a'] = 2*(Mgam/Sa+dsa)*Mgam*(1.-parmDict[phfx+'Mustrain:mx'])**2/ateln2       
    20962141    else:       #generalized - P.W. Stephens model
    20972142        pwrs = calcControls[phfx+'MuPwrs']
    20982143        const = 0.018*refl[4]**2*tanth
     2144        sum = 0
    20992145        for i,pwr in enumerate(pwrs):
    2100             gamDict[phfx+'Mustrain:'+str(i)] = const*refl[0]**pwr[0]*refl[1]**pwr[1]*refl[2]**pwr[2]
    2101     return gamDict
     2146            term = refl[0]**pwr[0]*refl[1]**pwr[1]*refl[2]**pwr[2]
     2147            sum += parmDict[phfx+'Mustrain:'+str(i)]*term
     2148            gamDict[phfx+'Mustrain:'+str(i)] = const*term*parmDict[phfx+'Mustrain:mx']
     2149            sigDict[phfx+'Mustrain:'+str(i)] = \
     2150                2.*const*term*(1.-parmDict[phfx+'Mustrain:mx'])**2/ateln2
     2151        Mgam = 0.018*refl[4]**2*tand(refl[5]/2.)*sum
     2152        for i in range(len(pwrs)):
     2153            sigDict[phfx+'Mustrain:'+str(i)] *= Mgam
     2154    gamDict[phfx+'Mustrain:mx'] = Mgam
     2155    sigDict[phfx+'Mustrain:mx'] = -2.*Mgam**2*(1.-parmDict[phfx+'Mustrain:mx'])/ateln2
     2156    return sigDict,gamDict
    21022157       
    21032158def GetReflPos(refl,wave,G,hfx,calcControls,parmDict):
     
    22762331def getPowderProfile(parmDict,x,varylist,Histogram,Phases,calcControls,pawleyLookup):
    22772332   
    2278     def GetReflSIgGam(refl,wave,G,GB,hfx,phfx,calcControls,parmDict):
     2333    def GetReflSigGam(refl,wave,G,GB,hfx,phfx,calcControls,parmDict):
    22792334        U = parmDict[hfx+'U']
    22802335        V = parmDict[hfx+'V']
     
    22832338        Y = parmDict[hfx+'Y']
    22842339        tanPos = tand(refl[5]/2.0)
    2285         sig = U*tanPos**2+V*tanPos+W        #save peak sigma
     2340        Ssig,Sgam = GetSampleSigGam(refl,wave,G,GB,phfx,calcControls,parmDict)
     2341        sig = U*tanPos**2+V*tanPos+W+Ssig     #save peak sigma
    22862342        sig = max(0.001,sig)
    2287         gam = X/cosd(refl[5]/2.0)+Y*tanPos+GetSampleGam(refl,wave,G,GB,phfx,calcControls,parmDict) #save peak gamma
     2343        gam = X/cosd(refl[5]/2.0)+Y*tanPos+Sgam    #save peak gamma
    22882344        gam = max(0.001,gam)
    22892345        return sig,gam
     
    23282384                Lorenz = 1./(2.*sind(refl[5]/2.)**2*cosd(refl[5]/2.))           #Lorentz correction
    23292385                refl[5] += GetHStrainShift(refl,SGData,phfx,parmDict)               #apply hydrostatic strain shift
    2330                 refl[6:8] = GetReflSIgGam(refl,wave,G,GB,hfx,phfx,calcControls,parmDict)    #peak sig & gam
     2386                refl[6:8] = GetReflSigGam(refl,wave,G,GB,hfx,phfx,calcControls,parmDict)    #peak sig & gam
    23312387                GetIntensityCorr(refl,G,g,pfx,phfx,hfx,SGData,calcControls,parmDict)    #puts corrections in refl[13]
    23322388                refl[13] *= Vst*Lorenz
     
    25492605                        if Ka2:
    25502606                            depDerivDict[name][iBeg2:iFin2] += dDijDict[name]*dervDict2['pos']
    2551                 gamDict = GetSampleGamDerv(refl,wave,G,GB,phfx,calcControls,parmDict)
     2607                sigDict,gamDict = GetSampleSigGamDerv(refl,wave,G,GB,phfx,calcControls,parmDict)
    25522608                for name in gamDict:
    25532609                    if name in varylist:
     
    25592615                        if Ka2:
    25602616                            depDerivDict[name][iBeg2:iFin2] += gamDict[name]*dervDict2['gam']
     2617                for name in sigDict:
     2618                    if name in varylist:
     2619                        dMdv[varylist.index(name)][iBeg:iFin] += sigDict[name]*dervDict['sig']
     2620                        if Ka2:
     2621                            dMdv[varylist.index(name)][iBeg2:iFin2] += sigDict[name]*dervDict2['sig']
     2622                    elif name in dependentVars:
     2623                        depDerivDict[name][iBeg:iFin] += sigDict[name]*dervDict['sig']
     2624                        if Ka2:
     2625                            depDerivDict[name][iBeg2:iFin2] += sigDict[name]*dervDict2['sig']
    25612626                                               
    25622627            elif 'T' in calcControls[hfx+'histType']:
     
    26262691    return Vec,Hess
    26272692
    2628 #def HessRefine(values,HistoPhases,parmdict,varylist,calcControls,pawleyLookup,dlg):
    2629 #    parmdict.update(zip(varylist,values))
    2630 #    G2mv.Dict2Map(parmdict,varylist)
    2631 #    Histograms,Phases = HistoPhases
    2632 #    nvar = len(varylist)
    2633 #    Hess = np.empty(0)
    2634 #    histoList = Histograms.keys()
    2635 #    histoList.sort()
    2636 #    for histogram in histoList:
    2637 #        if 'PWDR' in histogram[:4]:
    2638 #            Histogram = Histograms[histogram]
    2639 #            hId = Histogram['hId']
    2640 #            hfx = ':%d:'%(hId)
    2641 #            Limits = calcControls[hfx+'Limits']
    2642 #            x,y,w,yc,yb,yd = Histogram['Data']
    2643 #            dy = y-yc
    2644 #            xB = np.searchsorted(x,Limits[0])
    2645 #            xF = np.searchsorted(x,Limits[1])
    2646 #            dMdvh = np.sqrt(w[xB:xF])*getPowderProfileDerv(parmdict,x[xB:xF],
    2647 #                varylist,Histogram,Phases,calcControls,pawleyLookup)
    2648 #            if dlg:
    2649 #                dlg.Update(Histogram['wRp'],newmsg='Hessian for histogram %d Rwp=%8.3f%s'%(hId,Histogram['wRp'],'%'))[0]
    2650 #            if len(Hess):
    2651 #                Vec += np.sum(dMdvh*np.sqrt(w[xB:xF])*dy[xB:xF],axis=1)
    2652 #                Hess += np.inner(dMdvh,dMdvh)
    2653 #            else:
    2654 #                Vec = np.sum(dMdvh*np.sqrt(w[xB:xF])*dy[xB:xF],axis=1)
    2655 #                Hess = np.inner(dMdvh,dMdvh)
    2656 #    return Vec,Hess
    2657 
    26582693def HessRefine(values,HistoPhases,parmdict,varylist,calcControls,pawleyLookup,dlg):
    26592694    parmdict.update(zip(varylist,values))
     
    26612696    Histograms,Phases = HistoPhases
    26622697    nvar = len(varylist)
    2663     HessSum = np.zeros((nvar,nvar))
    2664     VecSum = np.zeros(nvar)
     2698    Hess = np.empty(0)
    26652699    histoList = Histograms.keys()
    26662700    histoList.sort()
    2667     argList = []
    2668     MaxProcess = calcControls['max Hprocess']
    26692701    for histogram in histoList:
    26702702        if 'PWDR' in histogram[:4]:
    26712703            Histogram = Histograms[histogram]
    2672             argList.append([
    2673                 Histogram,parmdict,varylist,Phases,
    2674                 calcControls,pawleyLookup])
    2675     if MaxProcess > 1:
    2676         mpPool = mp.Pool(processes=MaxProcess)
    2677         results = mpPool.map(ComputePowderHessian,argList)
    2678         for Vec,Hess in results:
    2679             VecSum += Vec
    2680             HessSum += Hess
    2681     else:
    2682         for arg in argList:
    2683             Vec,Hess = ComputePowderHessian(arg)
    2684             VecSum += Vec
    2685             HessSum += Hess
    2686     return VecSum,HessSum
     2704            hId = Histogram['hId']
     2705            hfx = ':%d:'%(hId)
     2706            Limits = calcControls[hfx+'Limits']
     2707            x,y,w,yc,yb,yd = Histogram['Data']
     2708            dy = y-yc
     2709            xB = np.searchsorted(x,Limits[0])
     2710            xF = np.searchsorted(x,Limits[1])
     2711            dMdvh = np.sqrt(w[xB:xF])*getPowderProfileDerv(parmdict,x[xB:xF],
     2712                varylist,Histogram,Phases,calcControls,pawleyLookup)
     2713            if dlg:
     2714                dlg.Update(Histogram['wRp'],newmsg='Hessian for histogram %d Rwp=%8.3f%s'%(hId,Histogram['wRp'],'%'))[0]
     2715            if len(Hess):
     2716                Vec += np.sum(dMdvh*np.sqrt(w[xB:xF])*dy[xB:xF],axis=1)
     2717                Hess += np.inner(dMdvh,dMdvh)
     2718            else:
     2719                Vec = np.sum(dMdvh*np.sqrt(w[xB:xF])*dy[xB:xF],axis=1)
     2720                Hess = np.inner(dMdvh,dMdvh)
     2721    return Vec,Hess
    26872722
    26882723def ComputePowderProfile(args):
     
    26992734    return xB,xF,yc,yb,Histogram['Reflection Lists']
    27002735
    2701 #def errRefine(values,HistoPhases,parmdict,varylist,calcControls,pawleyLookup,dlg):       
    2702 #    parmdict.update(zip(varylist,values))
    2703 #    Values2Dict(parmdict, varylist, values)
    2704 #    G2mv.Dict2Map(parmdict,varylist)
    2705 #    Histograms,Phases = HistoPhases
    2706 #    M = np.empty(0)
    2707 #    sumwYo = 0
    2708 #    Nobs = 0
    2709 #    histoList = Histograms.keys()
    2710 #    histoList.sort()
    2711 #    for histogram in histoList:
    2712 #        if 'PWDR' in histogram[:4]:
    2713 #            Histogram = Histograms[histogram]
    2714 #            hId = Histogram['hId']
    2715 #            hfx = ':%d:'%(hId)
    2716 #            Limits = calcControls[hfx+'Limits']
    2717 #            x,y,w,yc,yb,yd = Histogram['Data']
    2718 #            yc *= 0.0                           #zero full calcd profiles
    2719 #            yb *= 0.0
    2720 #            yd *= 0.0
    2721 #            xB = np.searchsorted(x,Limits[0])
    2722 #            xF = np.searchsorted(x,Limits[1])
    2723 #            Histogram['Nobs'] = xF-xB
    2724 #            Nobs += Histogram['Nobs']
    2725 #            Histogram['sumwYo'] = np.sum(w[xB:xF]*y[xB:xF]**2)
    2726 #            sumwYo += Histogram['sumwYo']
    2727 #            yc[xB:xF],yb[xB:xF] = getPowderProfile(parmdict,x[xB:xF],
    2728 #                varylist,Histogram,Phases,calcControls,pawleyLookup)
    2729 #            yc[xB:xF] += yb[xB:xF]
    2730 #            yd[xB:xF] = y[xB:xF]-yc[xB:xF]
    2731 #            Histogram['sumwYd'] = np.sum(np.sqrt(w[xB:xF])*(yd[xB:xF]))
    2732 #            wdy = -np.sqrt(w[xB:xF])*(yd[xB:xF])
    2733 #            Histogram['wRp'] = min(100.,np.sqrt(np.sum(wdy**2)/Histogram['sumwYo'])*100.)
    2734 #            if dlg:
    2735 #                dlg.Update(Histogram['wRp'],newmsg='For histogram %d Rwp=%8.3f%s'%(hId,Histogram['wRp'],'%'))[0]
    2736 #            M = np.concatenate((M,wdy))
    2737 #    Histograms['sumwYo'] = sumwYo
    2738 #    Histograms['Nobs'] = Nobs
    2739 #    Rwp = min(100.,np.sqrt(np.sum(M**2)/sumwYo)*100.)
    2740 #    if dlg:
    2741 #        GoOn = dlg.Update(Rwp,newmsg='%s%8.3f%s'%('Powder profile Rwp =',Rwp,'%'))[0]
    2742 #        if not GoOn:
    2743 #            parmDict['saved values'] = values
    2744 #            raise Exception         #Abort!!
    2745 #    return M
    2746 #   
    27472736def errRefine(values,HistoPhases,parmdict,varylist,calcControls,pawleyLookup,dlg):       
    27482737    parmdict.update(zip(varylist,values))
     
    27552744    histoList = Histograms.keys()
    27562745    histoList.sort()
    2757     argList = []
    2758     MaxProcess = calcControls['max Hprocess']
    27592746    for histogram in histoList:
    27602747        if 'PWDR' in histogram[:4]:
    27612748            Histogram = Histograms[histogram]
    2762             argList.append(
    2763                 [Histogram,parmdict,varylist,Phases,calcControls,pawleyLookup]
    2764                 )
    2765     if MaxProcess > 1:
    2766         mpPool = mp.Pool(processes=MaxProcess)
    2767         results = mpPool.map(ComputePowderProfile,argList)
    2768         for arg,res in zip(argList,results):
    2769             xB,xF,ycSect,ybSect,RL = res
    2770             Histogram = arg[0]
    2771             Histogram['Reflection Lists'] = RL
     2749            hId = Histogram['hId']
     2750            hfx = ':%d:'%(hId)
     2751            Limits = calcControls[hfx+'Limits']
    27722752            x,y,w,yc,yb,yd = Histogram['Data']
    27732753            yc *= 0.0                           #zero full calcd profiles
    27742754            yb *= 0.0
    27752755            yd *= 0.0
     2756            xB = np.searchsorted(x,Limits[0])
     2757            xF = np.searchsorted(x,Limits[1])
    27762758            Histogram['Nobs'] = xF-xB
    27772759            Nobs += Histogram['Nobs']
    27782760            Histogram['sumwYo'] = np.sum(w[xB:xF]*y[xB:xF]**2)
    27792761            sumwYo += Histogram['sumwYo']
    2780            
    2781             yc[xB:xF] = ycSect
    2782             yb[xB:xF] = ybSect
    2783             yc[xB:xF] += yb[xB:xF]
    2784             yd[xB:xF] = y[xB:xF]-yc[xB:xF]
    2785             Histogram['sumwYd'] = np.sum(np.sqrt(w[xB:xF])*(yd[xB:xF]))
    2786             wdy = -np.sqrt(w[xB:xF])*(yd[xB:xF])
    2787             Histogram['wRp'] = min(100.,np.sqrt(np.sum(wdy**2)/Histogram['sumwYo'])*100.)
    2788             M = np.concatenate((M,wdy))
    2789     else:
    2790         for arg in argList:
    2791             xB,xF,ycSect,ybSect,RL = ComputePowderProfile(arg)
    2792             Histogram = arg[0]
    2793             hId = Histogram['hId']
    2794             x,y,w,yc,yb,yd = Histogram['Data']
    2795             yc *= 0.0                           #zero full calcd profiles
    2796             yb *= 0.0
    2797             yd *= 0.0
    2798             Histogram['Nobs'] = xF-xB
    2799             Nobs += Histogram['Nobs']
    2800             Histogram['sumwYo'] = np.sum(w[xB:xF]*y[xB:xF]**2)
    2801             sumwYo += Histogram['sumwYo']
    2802            
    2803             yc[xB:xF] = ycSect
    2804             yb[xB:xF] = ybSect
     2762            yc[xB:xF],yb[xB:xF] = getPowderProfile(parmdict,x[xB:xF],
     2763                varylist,Histogram,Phases,calcControls,pawleyLookup)
    28052764            yc[xB:xF] += yb[xB:xF]
    28062765            yd[xB:xF] = y[xB:xF]-yc[xB:xF]
     
    28112770                dlg.Update(Histogram['wRp'],newmsg='For histogram %d Rwp=%8.3f%s'%(hId,Histogram['wRp'],'%'))[0]
    28122771            M = np.concatenate((M,wdy))
    2813 
    28142772    Histograms['sumwYo'] = sumwYo
    28152773    Histograms['Nobs'] = Nobs
     
    28202778            parmDict['saved values'] = values
    28212779            raise Exception         #Abort!!
    2822     return M       
    2823                    
     2780    return M
     2781                       
    28242782def Refine(GPXfile,dlg):
    28252783    import pytexture as ptx
     
    28722830    G2mv.Map2Dict(parmDict,varyList)
    28732831#    print G2mv.VarRemapShow(varyList)
    2874 
     2832    Rvals = {}
    28752833    while True:
    28762834        begin = time.time()
     
    28872845            result = G2mth.HessianLSQ(errRefine,values,Hess=HessRefine,ftol=Ftol,maxcyc=maxCyc,
    28882846                args=([Histograms,Phases],parmDict,varyList,calcControls,pawleyLookup,dlg))
    2889             ncyc = result[2]['num cyc']+1                           
     2847            ncyc = result[2]['num cyc']+1
     2848            Rvals['lamMax'] = result[2]['lamMax']                           
    28902849        else:           #'numeric'
    28912850            result = so.leastsq(errRefine,values,full_output=True,ftol=Ftol,epsfcn=1.e-8,factor=Factor,
     
    28952854#        for item in table: print item,table[item]               #useful debug - are things shifting?
    28962855        runtime = time.time()-begin
    2897         chisq = np.sum(result[2]['fvec']**2)
     2856        Rvals['chisq'] = np.sum(result[2]['fvec']**2)
    28982857        Values2Dict(parmDict, varyList, result[0])
    28992858        G2mv.Dict2Map(parmDict,varyList)
    29002859       
    2901         Rwp = np.sqrt(chisq/Histograms['sumwYo'])*100.      #to %
    2902         GOF = chisq/(Histograms['Nobs']-len(varyList))
     2860        Rvals['Nobs'] = Histograms['Nobs']
     2861        Rvals['Rwp'] = np.sqrt(Rvals['chisq']/Histograms['sumwYo'])*100.      #to %
     2862        Rvals['GOF'] = Rvals['chisq']/(Histograms['Nobs']-len(varyList))
    29032863        print '\n Refinement results:'
    29042864        print 135*'-'
    29052865        print ' Number of function calls:',result[2]['nfev'],' Number of observations: ',Histograms['Nobs'],' Number of parameters: ',len(varyList)
    29062866        print ' Refinement time = %8.3fs, %8.3fs/cycle, for %d cycles'%(runtime,runtime/ncyc,ncyc)
    2907         print ' wRp = %7.2f%%, chi**2 = %12.6g, reduced chi**2 = %6.2f'%(Rwp,chisq,GOF)
     2867        print ' wRp = %7.2f%%, chi**2 = %12.6g, reduced chi**2 = %6.2f'%(Rvals['Rwp'],Rvals['chisq'],Rvals['GOF'])
    29082868        try:
    2909             covMatrix = result[1]*GOF
     2869            covMatrix = result[1]*Rvals['GOF']
    29102870            sig = np.sqrt(np.diag(covMatrix))
    29112871            if np.any(np.isnan(sig)):
     
    29392899    newCellDict = GetNewCellParms(parmDict,varyList)
    29402900    newAtomDict = ApplyXYZshifts(parmDict,varyList)
    2941     covData = {'variables':result[0],'varyList':varyList,'sig':sig,
     2901    covData = {'variables':result[0],'varyList':varyList,'sig':sig,'Rvals':Rvals,
    29422902        'covMatrix':covMatrix,'title':GPXfile,'newAtomDict':newAtomDict,'newCellDict':newCellDict}
    29432903    # add the uncertainties into the esd dictionary (sigDict)
     
    29642924
    29652925    if dlg:
    2966         return Rwp
     2926        return Rvals['Rwp']
    29672927
    29682928def SeqRefine(GPXfile,dlg):
     
    30413001        G2mv.GenerateConstraints(groups,parmlist,varyList,constrDict,constrFlag,fixedList)
    30423002        G2mv.Map2Dict(parmDict,varyList)
    3043    
     3003        Rvals = {}
    30443004        while True:
    30453005            begin = time.time()
     
    30663026
    30673027            runtime = time.time()-begin
    3068             chisq = np.sum(result[2]['fvec']**2)
     3028            Rvals['chisq'] = np.sum(result[2]['fvec']**2)
    30693029            Values2Dict(parmDict, varyList, result[0])
    30703030            G2mv.Dict2Map(parmDict,varyList)
    30713031           
    3072             Rwp = np.sqrt(chisq/Histo['sumwYo'])*100.      #to %
    3073             GOF = chisq/(Histo['Nobs']-len(varyList))
     3032            Rvals['Rwp'] = np.sqrt(Rvals['chisq']/Histo['sumwYo'])*100.      #to %
     3033            Rvals['GOF'] = Rvals['Rwp']/(Histo['Nobs']-len(varyList))
     3034            Rvals['Nobs'] = Histo['Nobs']
    30743035            print '\n Refinement results for histogram: v'+histogram
    30753036            print 135*'-'
    30763037            print ' Number of function calls:',result[2]['nfev'],' Number of observations: ',Histo['Nobs'],' Number of parameters: ',len(varyList)
    30773038            print ' Refinement time = %8.3fs, %8.3fs/cycle, for %d cycles'%(runtime,runtime/ncyc,ncyc)
    3078             print ' wRp = %7.2f%%, chi**2 = %12.6g, reduced chi**2 = %6.2f'%(Rwp,chisq,GOF)
     3039            print ' wRp = %7.2f%%, chi**2 = %12.6g, reduced chi**2 = %6.2f'%(Rvals['Rwp'],Rvals['chisq'],Rvals['GOF'])
    30793040            try:
    30803041                covMatrix = result[1]*GOF
     
    31023063        newCellDict = GetNewCellParms(parmDict,varyList)
    31033064        newAtomDict = ApplyXYZshifts(parmDict,varyList)
    3104         covData = {'variables':result[0],'varyList':varyList,'sig':sig,
     3065        covData = {'variables':result[0],'varyList':varyList,'sig':sig,'Rvals':Rvals,
    31053066            'covMatrix':covMatrix,'title':histogram,'newAtomDict':newAtomDict,'newCellDict':newCellDict}
    31063067        SetHistogramPhaseData(parmDict,sigDict,Phases,Histo,ifPrint)
Note: See TracChangeset for help on using the changeset viewer.