Changeset 1459


Ignore:
Timestamp:
Aug 8, 2014 2:59:09 PM (7 years ago)
Author:
vondreele
Message:

add instrument parameters (flight path & detector 2-theta) needed for TOF
rework reflection list for TOF
change default diff curve & reflection marker offsets
change weighting to instrument constants calibration to be 1/esd2 from peak fit positions - works a lot better
1st shot at TOF Rietveld refinement with derivatives - need to be checked for correctness (some are wrong)

Location:
trunk
Files:
7 edited

Legend:

Unmodified
Added
Removed
  • trunk/GSASII.py

    r1453 r1459  
    979979                return [G2IO.makeInstDict(names,data,codes),{}]
    980980            elif 'T' in DataType:
    981                 names = ['Type','2-theta','difC','difA','difB','Zero','alpha','beta-0','beta-1',
    982                     'beta-q','sig-0','sig-1','sig-q','X','Y','Azimuth']
    983                 codes = [0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0]
     981                names = ['Type','fltPath','2-theta','difC','difA', 'difB','Zero','alpha','beta-0','beta-1',
     982                    'beta-q','sig-0','sig-1','sig-q','X', 'Y','Azimuth',]
     983                codes = [0,0,0,0,0, 0,0,0,0,0, 0,0,0,0,0, 0,0,]
    984984                azm = 0.
    985985                if 'INS  1DETAZM' in Iparm:
    986986                    azm = float(Iparm['INS  1DETAZM'])
     987                s = Iparm['INS   FPATH1'].split()
     988                fltPath0 = G2IO.sfloat(s[0])
    987989                s = Iparm['INS  1BNKPAR'].split()
     990                fltPath1 = G2IO.sfloat(s[0])
     991                data.extend([fltPath0+fltPath1,])               #Flight path source-sample-detector
    988992                data.extend([G2IO.sfloat(s[1]),])               #2-theta for bank
    989993                s = Iparm['INS  1 ICONS'].split()
     
    18891893        self.Offset = [0.0,0.0]
    18901894        self.delOffset = .02
    1891         self.refOffset = -100.0
     1895        self.refOffset = -1.0
    18921896        self.refDelt = .01
    18931897        self.Weight = False
     
    27912795                            if name2 == 'Peak List':
    27922796                                peaks = self.PatternTree.GetItemPyData(item2)['peaks']
    2793                                 file.write("%s \n" % (name+' Peak List'))               
     2797                                file.write("%s \n" % (name+' Peak List'))
     2798                                if len(peaks[0]) == 8:
     2799                                    file.write('%10s %12s %10s %10s\n'%('pos','int','sig','gam'))
     2800                                else:
     2801                                    file.write('%10s %12s %10s %10s %10s %10s\n'%('pos','int','alp','bet','sig','gam'))                                   
    27942802                                for peak in peaks:
    2795                                     file.write("%10.5f %12.2f %10.3f %10.3f \n" % \
    2796                                         (peak[0],peak[2],peak[4],peak[6]))
     2803                                    if len(peak) == 8:  #CW
     2804                                        file.write("%10.5f %12.2f %10.3f %10.3f \n" % \
     2805                                            (peak[0],peak[2],peak[4],peak[6]))
     2806                                    else:               #TOF - more cols
     2807                                        file.write("%10.5f %12.2f %10.3f %10.3f %10.3f %10.3f\n" % \
     2808                                            (peak[0],peak[2],peak[4],peak[6],peak[8],peak[10]))
    27972809                            item2, cookie2 = self.PatternTree.GetNextChild(item, cookie2)                           
    27982810                    item, cookie = self.PatternTree.GetNextChild(self.root, cookie)                           
  • trunk/GSASIIIO.py

    r1446 r1459  
    16941694    defaultIparm_lbl.append('10m TOF backscattering bank')
    16951695    defaultIparms.append({
     1696        'INS   FPATH1':'      9.00',
    16961697        'INS   HTYPE ':'PNT',
    16971698        'INS  1 ICONS':'   5000.00      0.00      0.00',
     
    17031704    defaultIparm_lbl.append('10m TOF 90deg bank')
    17041705    defaultIparms.append({
     1706        'INS   FPATH1':'      9.00',
    17051707        'INS   HTYPE ':'PNT',
    17061708        'INS  1 ICONS':'   3500.00      0.00      0.00',
     
    17121714    defaultIparm_lbl.append('63m POWGEN 90deg bank')
    17131715    defaultIparms.append({
     1716        'INS   FPATH1':'     60.00',
    17141717        'INS   HTYPE ':'PNT',
    17151718        'INS  1 ICONS':'  22585.80      0.00      0.00',
    1716         'INS  1BNKPAR':'    1.0000    90.000',       
     1719        'INS  1BNKPAR':'     3.169    90.000',       
    17171720        'INS  1PRCF1 ':'    1    8   0.01000',
    17181721        'INS  1PRCF11':'   0.000000E+00   1.000000E+00   3.000000E-02   4.000000E-03',
  • trunk/GSASIIplot.py

    r1451 r1459  
    874874                if G2frame.SqrtPlot:
    875875                    G2frame.delOffset = .002
    876                     G2frame.refOffset = -10.0
     876                    G2frame.refOffset = -1.0
    877877                    G2frame.refDelt = .001
    878878                else:
    879879                    G2frame.delOffset = .02
    880                     G2frame.refOffset = -100.0
     880                    G2frame.refOffset = -1.0
    881881                    G2frame.refDelt = .01
    882882            newPlot = True
     
    18861886    else:
    18871887        Plot.set_ylabel(r'$\mathsf{\Delta}T/T$',fontsize=14)
    1888     colors=['b','g','r','c','m','k']
    18891888    for ixy,xy in enumerate(XY):
    18901889        X,Y = xy
     
    18991898            Plot.errorbar(X,Y,ecolor='k',yerr=E)
    19001899        Plot.plot(X,Y,'kx',picker=3)
     1900        Plot.axhline(0.,color='r',linestyle='--')
    19011901    if not newPlot:
    19021902        Page.toolbar.push_current()
  • trunk/GSASIIpwd.py

    r1443 r1459  
    561561            bakInt = si.interp1d(bakPos,bakVals,'linear')
    562562            yb = bakInt(xdata)
    563     if 'difC' in parmDict:
     563    if pfx+'difC' in parmDict:
    564564        ff = 1.
    565565    else:       
     
    602602    return yb
    603603   
    604 def getBackgroundDerv(pfx,parmDict,bakType,xdata):
     604def getBackgroundDerv(hfx,parmDict,bakType,xdata):
    605605    'needs a doc string'
    606606    nBak = 0
    607607    while True:
    608         key = pfx+'Back:'+str(nBak)
     608        key = hfx+'Back:'+str(nBak)
    609609        if key in parmDict:
    610610            nBak += 1
     
    612612            break
    613613    dydb = np.zeros(shape=(nBak,len(xdata)))
    614     dyddb = np.zeros(shape=(3*parmDict[pfx+'nDebye'],len(xdata)))
    615     dydpk = np.zeros(shape=(4*parmDict[pfx+'nPeaks'],len(xdata)))
     614    dyddb = np.zeros(shape=(3*parmDict[hfx+'nDebye'],len(xdata)))
     615    dydpk = np.zeros(shape=(4*parmDict[hfx+'nPeaks'],len(xdata)))
    616616    cw = np.diff(xdata)
    617617    cw = np.append(cw,cw[-1])
     
    650650                        np.where(xdata<bakPos[i+1],(bakPos[i+1]-xdata)/(bakPos[i+1]-bakPos[i]),0.),
    651651                        np.where(xdata>bakPos[i-1],(xdata-bakPos[i-1])/(bakPos[i]-bakPos[i-1]),0.))
    652     if 'difC' in parmDict:
     652    if hfx+'difC' in parmDict:
    653653        ff = 1.
    654654    else:
    655655        try:
    656             wave = parmDict[pfx+'Lam']
     656            wave = parmDict[hfx+'Lam']
    657657        except KeyError:
    658             wave = parmDict[pfx+'Lam1']
     658            wave = parmDict[hfx+'Lam1']
    659659        q = 4.0*np.pi*npsind(xdata/2.0)/wave
    660660        SQ = (q/(4*np.pi))**2
     
    664664    while True:
    665665        try:
    666             dbA = parmDict[pfx+'DebyeA:'+str(iD)]
    667             dbR = parmDict[pfx+'DebyeR:'+str(iD)]
    668             dbU = parmDict[pfx+'DebyeU:'+str(iD)]
     666            if hfx+'difC' in parmDict:
     667                q = 2*np.pi*parmDict[hfx+'difC']/xdata
     668            dbA = parmDict[hfx+'DebyeA:'+str(iD)]
     669            dbR = parmDict[hfx+'DebyeR:'+str(iD)]
     670            dbU = parmDict[hfx+'DebyeU:'+str(iD)]
    669671            sqr = np.sin(q*dbR)/(q*dbR)
    670672            cqr = np.cos(q*dbR)
     
    679681    while True:
    680682        try:
    681             pkP = parmDict[pfx+'BkPkpos;'+str(iD)]
    682             pkI = parmDict[pfx+'BkPkint;'+str(iD)]
    683             pkS = parmDict[pfx+'BkPksig;'+str(iD)]
    684             pkG = parmDict[pfx+'BkPkgam;'+str(iD)]
     683            pkP = parmDict[hfx+'BkPkpos;'+str(iD)]
     684            pkI = parmDict[hfx+'BkPkint;'+str(iD)]
     685            pkS = parmDict[hfx+'BkPksig;'+str(iD)]
     686            pkG = parmDict[hfx+'BkPkgam;'+str(iD)]
    685687            shl = 0.002
    686688            Wd,fmin,fmax = getWidthsCW(pkP,pkS,pkG,shl)
     
    11721174    peakDsp = []
    11731175    peakWt = []
    1174     for peak in IndexPeaks:
     1176    for peak,sig in zip(IndexPeaks[0],IndexPeaks[1]):
    11751177        if peak[2] and peak[3]:
    11761178            peakPos.append(peak[0])
    11771179            peakDsp.append(peak[8])
    1178             peakWt.append(1/peak[1])
     1180            peakWt.append(1/sig**2)
    11791181    peakPos = np.array(peakPos)
    11801182    peakDsp = np.array(peakDsp)
  • trunk/GSASIIpwdGUI.py

    r1453 r1459  
    8787        'SlitLen':0.0,                          #Slit length - in Q(A-1)
    8888        }
    89 def SetupSampleLabels(histName,dataType):
     89def SetupSampleLabels(histName,dataType,histType):
    9090    '''Setup a list of labels and number formatting for use in
    9191    labeling sample parameters.
     
    9595    parms = []
    9696    parms.append(['Scale','Histogram scale factor: ',[10,4]])
    97     parms.append(['Gonio. radius','Goniometer radius (mm): ',[10,3]])
     97    if 'C' in histType:
     98        parms.append(['Gonio. radius','Goniometer radius (mm): ',[10,3]])
    9899    if 'PWDR' in histName:
    99100        if dataType == 'Debye-Scherrer':
     
    101102                ['DisplaceY',u'Sample Y displ. || to beam (\xb5m): ',[10,3]],
    102103                ['Absorption',u'Sample absorption (\xb5\xb7r): ',[10,4]],]
     104            if 'T' in histType:
     105                parms[-1] = ['Absorption',u'Sample absorption (\xb5\xb7r/l): ',[10,4]]
    103106        elif dataType == 'Bragg-Brentano':
    104107            parms += [['Shift',u'Sample displacement(\xb5m): ',[10,4]],
     
    927930            if key in ['Type','U','V','W','X','Y','SH/L','I(L2)/I(L1)','alpha',
    928931                'beta-0','beta-1','beta-q','sig-0','sig-1','sig-q','Polariz.',
    929                 'Lam','Azimuth','2-theta','difC','difA','difB','Zero','Lam1','Lam2']:
     932                'Lam','Azimuth','2-theta','fltPath','difC','difA','difB','Zero','Lam1','Lam2']:
    930933                good.append(key)
    931934        return good
     
    964967            G2frame.ErrorDialog('Can not calibrate','Index Peak List not indexed')
    965968            return           
    966         G2pwd.DoCalibInst(IndexPeaks[0],data)
     969        G2pwd.DoCalibInst(IndexPeaks,data)
    967970        UpdateInstrumentGrid(G2frame,data)
    968971        XY = []
     
    12551258            else:                                   #time of flight (neutrons)
    12561259                subSizer = wx.BoxSizer(wx.HORIZONTAL)
     1260                subSizer.Add(wx.StaticText(G2frame.dataDisplay,-1,' Fligth path: '),0,WACV)
     1261                txt = '%8.3f'%(insVal['fltPath'])
     1262                subSizer.Add(wx.StaticText(G2frame.dataDisplay,-1,txt.strip()),0,WACV)
     1263                labelLst.append('flight path')
     1264                elemKeysLst.append(['fltpath',1])
     1265                dspLst.append([10,2])
     1266                refFlgElem.append(None)                   
    12571267                subSizer.Add(wx.StaticText(G2frame.dataDisplay,-1,'  2-theta: '),0,WACV)
    12581268                txt = '%7.2f'%(insVal['2-theta'])
     
    12731283                    subSizer.Add(wx.StaticText(G2frame.dataDisplay,-1,'  alpha, beta: fixed by table'),0,WACV)
    12741284                else:
    1275                     Items = ['difC','difA','difB','Zero','alpha','beta-0','beta-1','beta-q','sig-0','sig-1','sig-q','X','Y']
     1285                    Items = ['difC','difA','difB','Zero','alpha','beta-0','beta-1','beta-q','sig-0','sig-1','sig-q','','X','Y']
     1286                mainSizer.Add((5,5),0)
    12761287                mainSizer.Add(subSizer)
     1288                mainSizer.Add((5,5),0)
    12771289                for item in Items:
     1290                    if item == '':
     1291                        instSizer.Add((5,5),0)
     1292                        instSizer.Add((5,5),0)
     1293                        instSizer.Add((5,5),0)
     1294                        continue
    12781295                    nDig = (10,3)
    12791296                    fmt = '%10.3f'
     
    12921309                    refFlgElem.append([item,2])
    12931310                    instSizer.Add(RefineBox(item),0,WACV)
    1294 #                    if not ifHisto and item in ['difC','difA','difB','Zero',]:
    1295 #                        refFlgElem.append(None)
    1296 #                        instSizer.Add((5,5),0)
    1297 #                    else:
    1298 #                        refFlgElem.append([item,2])
    1299 #                        instSizer.Add(RefineBox(item),0,WACV)
    13001311        elif 'S' in insVal['Type']:                       #single crystal data
    13011312            if 'C' in insVal['Type']:               #constant wavelength
     
    15271538        # Assemble a list of item labels
    15281539        TextTable = {key:label for key,label,dig in
    1529                      SetupSampleLabels(hst,data.get('Type'))
     1540                     SetupSampleLabels(hst,data.get('Type'),Inst['Type'][0])
    15301541                     }
    15311542        # get flexible labels
     
    16481659    #reload(G2gd)
    16491660    ######################################################################
     1661    Inst = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(
     1662            G2frame,G2frame.PatternId, 'Instrument Parameters'))[0]
    16501663    histName = G2frame.PatternTree.GetItemText(G2frame.PatternId)
    16511664    if G2frame.dataDisplay:
     
    16981711#patch end
    16991712    labelLst,elemKeysLst,dspLst,refFlgElem = [],[],[],[]
    1700     parms = SetupSampleLabels(histName,data.get('Type'))
     1713    parms = SetupSampleLabels(histName,data.get('Type'),Inst['Type'][0])
    17011714    mainSizer = wx.BoxSizer(wx.VERTICAL)
    17021715    topSizer = wx.BoxSizer(wx.HORIZONTAL)
     
    17241737        nameSizer.Add(wx.StaticText(G2frame.dataDisplay,wx.ID_ANY,' Diffractometer type: '),
    17251738                    0,WACV)
    1726         choices = ['Debye-Scherrer','Bragg-Brentano',]
     1739        if 'T' in Inst['Type'][0]:
     1740            choices = ['Debye-Scherrer',]
     1741        else:
     1742            choices = ['Debye-Scherrer','Bragg-Brentano',]
    17271743        histoType = wx.ComboBox(G2frame.dataDisplay,wx.ID_ANY,value=data['Type'],choices=choices,
    17281744            style=wx.CB_READONLY|wx.CB_DROPDOWN)
     
    21232139        ibrav = bravaisSymb.index(controls[5])
    21242140        SGData = G2spc.SpcGroup(controls[13])[1]
    2125         dmin = G2indx.getDmin(peaks)-0.005
     2141        dmin = G2indx.getDmin(peaks[0])-0.005
    21262142        G2frame.HKL = G2pwd.getHKLpeak(dmin,SGData,A)
    21272143        G2indx.IndexPeaks(peaks[0],G2frame.HKL)
     
    24752491    if G2frame.dataDisplay:
    24762492        G2frame.dataFrame.Clear()
     2493    Inst = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Instrument Parameters'))[0]
    24772494    rowLabels = []
    24782495    if HKLF:
     
    24942511        except TypeError:
    24952512            refList = np.array([refl[:11] for refl in data[G2frame.RefList]])
    2496             I100 = refList.T[8]*np.array([refl[13] for refl in data[G2frame.RefList]])
     2513            I100 = refList.T[8]*np.array([refl[11] for refl in data[G2frame.RefList]])
    24972514        Imax = np.max(I100)
    24982515        if Imax:
    24992516            I100 *= 100.0/Imax
    2500         refs = np.vstack((refList.T[:11],I100)).T
     2517        if 'C' in Inst['Type'][0]:
     2518            refs = np.vstack((refList.T[:11],I100)).T
     2519        elif 'T' in Inst['Type'][0]:
     2520            refs = np.vstack((refList.T[:15],I100)).T
     2521           
    25012522    for i in range(len(refs)): rowLabels.append(str(i))
     2523    Types = 4*[wg.GRID_VALUE_LONG,]+4*[wg.GRID_VALUE_FLOAT+':10,4',]+ \
     2524        2*[wg.GRID_VALUE_FLOAT+':10,2',]+[wg.GRID_VALUE_FLOAT+':10,3',]+ \
     2525        [wg.GRID_VALUE_FLOAT+':10,3',]
    25022526    if HKLF:
    25032527        colLabels = ['H','K','L','mul','d','Fosq','sig','Fcsq','FoTsq','FcTsq','phase','Ext',]
    25042528    else:
    2505         colLabels = ['H','K','L','mul','d','pos','sig','gam','Fosq','Fcsq','phase','I100',]
    2506     Types = 4*[wg.GRID_VALUE_LONG,]+4*[wg.GRID_VALUE_FLOAT+':10,4',]+ \
    2507         2*[wg.GRID_VALUE_FLOAT+':10,2',]+[wg.GRID_VALUE_FLOAT+':10,3',]+ \
    2508         [wg.GRID_VALUE_FLOAT+':10,3',]
     2529        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',]
     2532        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',]
     2535           
    25092536    G2frame.PeakTable = G2gd.Table(refs,rowLabels=rowLabels,colLabels=colLabels,types=Types)
    25102537    G2frame.dataFrame.SetLabel('Reflection List for '+phaseName)
  • trunk/GSASIIstrIO.py

    r1415 r1459  
    239239            hId = Histogram['hId']
    240240            hfx = ':%d:'%(hId)
    241             keV = controlDict[hfx+'keV']
    242             for El in FFtables:
    243                 Orbs = G2el.GetXsectionCoeff(El.split('+')[0].split('-')[0])
    244                 FP,FPP,Mu = G2el.FPcalc(Orbs, keV)
    245                 FFtables[El][hfx+'FP'] = FP
    246                 FFtables[El][hfx+'FPP'] = FPP               
     241            if 'X' in controlDict[hfx+'histType']:
     242                keV = controlDict[hfx+'keV']
     243                for El in FFtables:
     244                    Orbs = G2el.GetXsectionCoeff(El.split('+')[0].split('-')[0])
     245                    FP,FPP,Mu = G2el.FPcalc(Orbs, keV)
     246                    FFtables[El][hfx+'FP'] = FP
     247                    FFtables[El][hfx+'FPP'] = FPP               
    247248           
    248249def GetPhaseNames(GPXfile):
     
    17181719                limits = Histogram['Limits'][1]
    17191720                inst = Histogram['Instrument Parameters'][0]
    1720                 Zero = inst['Zero'][1]
     1721                Zero = inst['Zero'][0]
    17211722                if 'C' in inst['Type'][1]:
    17221723                    try:
     
    17251726                        wave = inst['Lam1'][1]
    17261727                    dmin = wave/(2.0*sind(limits[1]/2.0))
     1728                elif 'T' in inst['Type'][0]:
     1729                    dmin = limits[0]/inst['difC'][1]
    17271730                pfx = str(pId)+':'+str(hId)+':'
    17281731                for item in ['Scale','Extinction']:
     
    18101813                            pos = 2.0*asind(wave/(2.0*d))+Zero
    18111814                            if limits[0] < pos < limits[1]:
    1812                                 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])
    18131816                                Uniq.append(uniq)
    18141817                                Phi.append(phi)
    1815                         else:
    1816                             raise ValueError
     1818                        elif 'T' in inst['Type'][0]:
     1819                            pos = inst['difC'][1]*d+inst['difA'][1]*d**2+inst['difB'][1]*d**3+Zero
     1820                            if limits[0] < pos < limits[1]:
     1821                                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                                Uniq.append(uniq)
     1824                                Phi.append(phi)
    18171825                    Histogram['Reflection Lists'][phase] = {'RefList':np.array(refList),'FF':{}}
    18181826            elif 'HKLF' in histogram:
     
    22182226            if len(Inst[item]) > 2 and Inst[item][2]:
    22192227                insVary.append(insName)
    2220 #        instDict[pfx+'X'] = max(instDict[pfx+'X'],0.001)
    2221 #        instDict[pfx+'Y'] = max(instDict[pfx+'Y'],0.001)
    2222         instDict[pfx+'SH/L'] = max(instDict[pfx+'SH/L'],0.0005)
     2228        if 'C' in dataType:
     2229            instDict[pfx+'SH/L'] = max(instDict[pfx+'SH/L'],0.0005)
    22232230        return dataType,instDict,insVary
    22242231       
     
    22812288    def PrintInstParms(Inst):
    22822289        print >>pFile,'\n Instrument Parameters:'
    2283         ptlbls = ' name  :'
    2284         ptstr =  ' value :'
    2285         varstr = ' refine:'
    22862290        insKeys = Inst.keys()
    22872291        insKeys.sort()
    2288         for item in insKeys:
    2289             if item not in ['Type','Source']:
    2290                 ptlbls += '%12s' % (item)
    2291                 ptstr += '%12.6f' % (Inst[item][1])
    2292                 if item in ['Lam1','Lam2','Azimuth']:
    2293                     varstr += 12*' '
    2294                 else:
    2295                     varstr += '%12s' % (str(bool(Inst[item][2])))
    2296         print >>pFile,ptlbls
    2297         print >>pFile,ptstr
    2298         print >>pFile,varstr
     2292        iBeg = 0
     2293        Ok = True
     2294        while Ok:
     2295            ptlbls = ' name  :'
     2296            ptstr =  ' value :'
     2297            varstr = ' refine:'
     2298            iFin = min(iBeg+9,len(insKeys))
     2299            for item in insKeys[iBeg:iFin]:
     2300                if item not in ['Type','Source']:
     2301                    ptlbls += '%12s' % (item)
     2302                    ptstr += '%12.6f' % (Inst[item][1])
     2303                    if item in ['Lam1','Lam2','Azimuth','fltPath','2-theta',]:
     2304                        varstr += 12*' '
     2305                    else:
     2306                        varstr += '%12s' % (str(bool(Inst[item][2])))
     2307            print >>pFile,ptlbls
     2308            print >>pFile,ptstr
     2309            print >>pFile,varstr
     2310            iBeg = iFin
     2311            if iBeg == len(insKeys):
     2312                Ok = False
     2313            else:
     2314                print >>pFile,'\n'
    22992315       
    23002316    def PrintSampleParms(Sample):
     
    23472363            Type,instDict,insVary = GetInstParms(hId,Inst)
    23482364            controlDict[pfx+'histType'] = Type
    2349             if pfx+'Lam1' in instDict:
    2350                 controlDict[pfx+'keV'] = 12.397639/instDict[pfx+'Lam1']
    2351             else:
    2352                 controlDict[pfx+'keV'] = 12.397639/instDict[pfx+'Lam']           
     2365            if 'XC' in Type:
     2366                if pfx+'Lam1' in instDict:
     2367                    controlDict[pfx+'keV'] = 12.397639/instDict[pfx+'Lam1']
     2368                else:
     2369                    controlDict[pfx+'keV'] = 12.397639/instDict[pfx+'Lam']           
    23532370            histDict.update(instDict)
    23542371            histVary += insVary
     
    23832400            Inst = Histogram['Instrument Parameters'][0]
    23842401            controlDict[pfx+'histType'] = Inst['Type'][0]
    2385             histDict[pfx+'Lam'] = Inst['Lam'][1]
    2386             controlDict[pfx+'keV'] = 12.397639/histDict[pfx+'Lam']                   
     2402            if 'X' in Inst['Type'][0]:
     2403                histDict[pfx+'Lam'] = Inst['Lam'][1]
     2404                controlDict[pfx+'keV'] = 12.397639/histDict[pfx+'Lam']                   
    23872405    return histVary,histDict,controlDict
    23882406   
     
    24982516       
    24992517    def PrintInstParmsSig(Inst,instSig):
    2500         ptlbls = ' names :'
    2501         ptstr =  ' value :'
    2502         sigstr = ' sig   :'
    25032518        refine = False
    25042519        insKeys = instSig.keys()
    25052520        insKeys.sort()
    2506         for name in insKeys:
    2507             if name not in  ['Type','Lam1','Lam2','Azimuth','Source']:
    2508                 ptlbls += '%12s' % (name)
    2509                 ptstr += '%12.6f' % (Inst[name][1])
    2510                 if instSig[name]:
    2511                     refine = True
    2512                     sigstr += '%12.6f' % (instSig[name])
    2513                 else:
    2514                     sigstr += 12*' '
    2515         if refine:
    2516             print >>pFile,'\n Instrument Parameters:'
    2517             print >>pFile,ptlbls
    2518             print >>pFile,ptstr
    2519             print >>pFile,sigstr
     2521        iBeg = 0
     2522        Ok = True
     2523        while Ok:
     2524            ptlbls = ' names :'
     2525            ptstr =  ' value :'
     2526            sigstr = ' sig   :'
     2527            iFin = min(iBeg+9,len(insKeys))
     2528            for name in insKeys[iBeg:iFin]:
     2529                if name not in  ['Type','Lam1','Lam2','Azimuth','Source','fltPath']:
     2530                    ptlbls += '%12s' % (name)
     2531                    ptstr += '%12.6f' % (Inst[name][1])
     2532                    if instSig[name]:
     2533                        refine = True
     2534                        sigstr += '%12.6f' % (instSig[name])
     2535                    else:
     2536                        sigstr += 12*' '
     2537            if refine:
     2538                print >>pFile,'\n Instrument Parameters:'
     2539                print >>pFile,ptlbls
     2540                print >>pFile,ptstr
     2541                print >>pFile,sigstr
     2542            iBeg = iFin
     2543            if iBeg == len(insKeys):
     2544                Ok = False
    25202545       
    25212546    def PrintSampleParmsSig(Sample,sampleSig):
  • trunk/GSASIIstrMath.py

    r1456 r1459  
    573573    for iref,refl in enumerate(refDict['RefList']):
    574574        if 'NT' in calcControls[hfx+'histType']:
    575             FP,FPP = G2el.BlenResCW(Tdata,BLtables,refl[12])
     575            FP,FPP = G2el.BlenResCW(Tdata,BLtables,refl[14])
    576576        fbs = np.array([0,0])
    577577        H = refl[:3]
     
    665665        SQfactor = 4.0*SQ*twopisq
    666666        if 'T' in calcControls[hfx+'histType']:
    667             FP,FPP = G2el.BlenResTOF(Tdata,BLtables,refl.T[12])
     667            if 'P' in calcControls[hfx+'histType']:
     668                FP,FPP = G2el.BlenResTOF(Tdata,BLtables,refl.T[14])
     669            else:
     670                FP,FPP = G2el.BlenResTOF(Tdata,BLtables,refl.T[12])
    668671            FP = np.repeat(FP.T,len(SGT),axis=0)
    669672            FPP = np.repeat(FPP.T,len(SGT),axis=0)
     
    710713    mSize = len(Mdata)
    711714    FF = np.zeros(len(Tdata))
    712     if 'N' in calcControls[hfx+'histType']:
     715    if 'NC' in calcControls[hfx+'histType']:
    713716        FP,FPP = G2el.BlenResCW(Tdata,BLtables,parmDict[hfx+'Lam'])
    714     else:
     717    elif 'X' in calcControls[hfx+'histType']:
    715718        FP = np.array([FFtables[El][hfx+'FP'] for El in Tdata])
    716719        FPP = np.array([FFtables[El][hfx+'FPP'] for El in Tdata])
     
    724727    dFdbab = np.zeros((nRef,2))
    725728    for iref,refl in enumerate(refDict['RefList']):
     729        if 'T' in calcControls[hfx+'histType']:
     730            FP,FPP = G2el.BlenResCW(Tdata,BLtables,refl.T[12])
    726731        H = np.array(refl[:3])
    727732        SQ = 1./(2.*refl[4])**2             # or (sin(theta)/lambda)**2
     
    731736        Tindx = np.array([refDict['FF']['El'].index(El) for El in Tdata])
    732737        FF = refDict['FF']['FF'][iref].T[Tindx]
    733 #        FF = [refDict['FF'][iref][El] for El in Tdata]         
    734738        Uniq = np.inner(H,SGMT)
    735739        Phi = np.inner(H,SGT)
     
    911915    'Spherical harmonics texture'
    912916    IFCoup = 'Bragg' in calcControls[hfx+'instType']
     917    if 'T' in calcControls[hfx+'histType']:
     918        tth = parmDict[hfx+'2-theta']
     919    else:
     920        tth = refl[5]
    913921    odfCor = 1.0
    914922    H = refl[:3]
     
    917925    Gangls = [parmDict[hfx+'Omega'],parmDict[hfx+'Chi'],parmDict[hfx+'Phi'],parmDict[hfx+'Azimuth']]
    918926    phi,beta = G2lat.CrsAng(H,cell,SGData)
    919     psi,gam,x,x = G2lat.SamAng(refl[5]/2.,Gangls,Sangls,IFCoup) #ignore 2 sets of angle derivs.
     927    psi,gam,x,x = G2lat.SamAng(tth/2.,Gangls,Sangls,IFCoup) #ignore 2 sets of angle derivs.
    920928    SHnames = G2lat.GenSHCoeff(SGData['SGLaue'],parmDict[pfx+'SHmodel'],parmDict[pfx+'SHorder'])
    921929    for item in SHnames:
     
    929937def SHTXcalDerv(refl,g,pfx,hfx,SGData,calcControls,parmDict):
    930938    'Spherical harmonics texture derivatives'
     939    if 'T' in calcControls[hfx+'histType']:
     940        tth = parmDict[hfx+'2-theta']
     941    else:
     942        tth = refl[5]
    931943    FORPI = 4.0*np.pi
    932944    IFCoup = 'Bragg' in calcControls[hfx+'instType']
     
    939951    Gangls = [parmDict[hfx+'Omega'],parmDict[hfx+'Chi'],parmDict[hfx+'Phi'],parmDict[hfx+'Azimuth']]
    940952    phi,beta = G2lat.CrsAng(H,cell,SGData)
    941     psi,gam,dPSdA,dGMdA = G2lat.SamAng(refl[5]/2.,Gangls,Sangls,IFCoup)
     953    psi,gam,dPSdA,dGMdA = G2lat.SamAng(tth/2.,Gangls,Sangls,IFCoup)
    942954    SHnames = G2lat.GenSHCoeff(SGData['SGLaue'],parmDict[pfx+'SHmodel'],parmDict[pfx+'SHorder'])
    943955    for item in SHnames:
     
    954966def SHPOcal(refl,g,phfx,hfx,SGData,calcControls,parmDict):
    955967    'spherical harmonics preferred orientation (cylindrical symmetry only)'
     968    if 'T' in calcControls[hfx+'histType']:
     969        tth = parmDict[hfx+'2-theta']
     970    else:
     971        tth = refl[5]
    956972    odfCor = 1.0
    957973    H = refl[:3]
     
    965981        IFCoup = False
    966982    phi,beta = G2lat.CrsAng(H,cell,SGData)
    967     psi,gam,x,x = G2lat.SamAng(refl[5]/2.,Gangls,Sangl,IFCoup) #ignore 2 sets of angle derivs.
     983    psi,gam,x,x = G2lat.SamAng(tth/2.,Gangls,Sangl,IFCoup) #ignore 2 sets of angle derivs.
    968984    SHnames = G2lat.GenSHCoeff(SGData['SGLaue'],'0',calcControls[phfx+'SHord'],False)
    969985    for item in SHnames:
     
    975991def SHPOcalDerv(refl,g,phfx,hfx,SGData,calcControls,parmDict):
    976992    'spherical harmonics preferred orientation derivatives (cylindrical symmetry only)'
     993    if 'T' in calcControls[hfx+'histType']:
     994        tth = parmDict[hfx+'2-theta']
     995    else:
     996        tth = refl[5]
    977997    FORPI = 12.5663706143592
    978998    odfCor = 1.0
     
    9881008        IFCoup = False
    9891009    phi,beta = G2lat.CrsAng(H,cell,SGData)
    990     psi,gam,x,x = G2lat.SamAng(refl[5]/2.,Gangls,Sangl,IFCoup) #ignore 2 sets of angle derivs.
     1010    psi,gam,x,x = G2lat.SamAng(tth/2.,Gangls,Sangl,IFCoup) #ignore 2 sets of angle derivs.
    9911011    SHnames = G2lat.GenSHCoeff(SGData['SGLaue'],'0',calcControls[phfx+'SHord'],False)
    9921012    for item in SHnames:
     
    10391059    'Needs a doc string'
    10401060    if 'Debye' in calcControls[hfx+'instType']:
    1041         return G2pwd.Absorb('Cylinder',parmDict[hfx+'Absorption'],refl[5],0,0)
     1061        if 'T' in calcControls[hfx+'histType']:
     1062            return G2pwd.Absorb('Cylinder',parmDict[hfx+'Absorption']*refl[14],parmDict[hfx+'2-theta'],0,0)
     1063        else:
     1064            return G2pwd.Absorb('Cylinder',parmDict[hfx+'Absorption'],refl[5],0,0)
    10421065    else:
    10431066        return G2pwd.SurfaceRough(parmDict[hfx+'SurfRoughA'],parmDict[hfx+'SurfRoughB'],refl[5])
     
    10461069    'Needs a doc string'
    10471070    if 'Debye' in calcControls[hfx+'instType']:
    1048         return G2pwd.AbsorbDerv('Cylinder',parmDict[hfx+'Absorption'],refl[5],0,0)
     1071        if 'T' in calcControls[hfx+'histType']:
     1072            return G2pwd.AbsorbDerv('Cylinder',parmDict[hfx+'Absorption']*refl[14],parmDict[hfx+'2-theta'],0,0)
     1073        else:
     1074            return G2pwd.AbsorbDerv('Cylinder',parmDict[hfx+'Absorption'],refl[5],0,0)
    10491075    else:
    10501076        return G2pwd.SurfaceRoughDerv(parmDict[hfx+'SurfRoughA'],parmDict[hfx+'SurfRoughB'],refl[5])
     
    10591085        Icorr *= SHTXcal(refl,g,pfx,hfx,SGData,calcControls,parmDict)
    10601086    Icorr *= GetAbsorb(refl,hfx,calcControls,parmDict)
    1061     refl[11] = Icorr       
     1087    return Icorr
    10621088   
    10631089def GetIntensityDerv(refl,uniq,G,g,pfx,phfx,hfx,SGData,calcControls,parmDict):
     
    10841110    return dIdsh,dIdsp,dIdPola,dIdPO,dFdODF,dFdSA,dFdAb
    10851111       
    1086 def GetSampleSigGam(refl,wave,G,GB,phfx,calcControls,parmDict):
     1112def GetSampleSigGam(refl,wave,G,GB,hfx,phfx,calcControls,parmDict):
    10871113    'Needs a doc string'
    1088     costh = cosd(refl[5]/2.)
    1089     #crystallite size
    1090     if calcControls[phfx+'SizeType'] == 'isotropic':
    1091         Sgam = 1.8*wave/(np.pi*parmDict[phfx+'Size;i']*costh)
    1092     elif calcControls[phfx+'SizeType'] == 'uniaxial':
    1093         H = np.array(refl[:3])
    1094         P = np.array(calcControls[phfx+'SizeAxis'])
    1095         cosP,sinP = G2lat.CosSinAngle(H,P,G)
    1096         Sgam = (1.8*wave/np.pi)/(parmDict[phfx+'Size;i']*parmDict[phfx+'Size;a']*costh)
    1097         Sgam *= np.sqrt((sinP*parmDict[phfx+'Size;a'])**2+(cosP*parmDict[phfx+'Size;i'])**2)
    1098     else:           #ellipsoidal crystallites
    1099         Sij =[parmDict[phfx+'Size:%d'%(i)] for i in range(6)]
    1100         H = np.array(refl[:3])
    1101         lenR = G2pwd.ellipseSize(H,Sij,GB)
    1102         Sgam = 1.8*wave/(np.pi*costh*lenR)
    1103     #microstrain               
    1104     if calcControls[phfx+'MustrainType'] == 'isotropic':
    1105         Mgam = 0.018*parmDict[phfx+'Mustrain;i']*tand(refl[5]/2.)/np.pi
    1106     elif calcControls[phfx+'MustrainType'] == 'uniaxial':
    1107         H = np.array(refl[:3])
    1108         P = np.array(calcControls[phfx+'MustrainAxis'])
    1109         cosP,sinP = G2lat.CosSinAngle(H,P,G)
    1110         Si = parmDict[phfx+'Mustrain;i']
    1111         Sa = parmDict[phfx+'Mustrain;a']
    1112         Mgam = 0.018*Si*Sa*tand(refl[5]/2.)/(np.pi*np.sqrt((Si*cosP)**2+(Sa*sinP)**2))
    1113     else:       #generalized - P.W. Stephens model
    1114         pwrs = calcControls[phfx+'MuPwrs']
    1115         sum = 0
    1116         for i,pwr in enumerate(pwrs):
    1117             sum += parmDict[phfx+'Mustrain:'+str(i)]*refl[0]**pwr[0]*refl[1]**pwr[1]*refl[2]**pwr[2]
    1118         Mgam = 0.018*refl[4]**2*tand(refl[5]/2.)*sum
     1114    if 'C' in calcControls[hfx+'histType']:
     1115        costh = cosd(refl[5]/2.)
     1116        #crystallite size
     1117        if calcControls[phfx+'SizeType'] == 'isotropic':
     1118            Sgam = 1.8*wave/(np.pi*parmDict[phfx+'Size;i']*costh)
     1119        elif calcControls[phfx+'SizeType'] == 'uniaxial':
     1120            H = np.array(refl[:3])
     1121            P = np.array(calcControls[phfx+'SizeAxis'])
     1122            cosP,sinP = G2lat.CosSinAngle(H,P,G)
     1123            Sgam = (1.8*wave/np.pi)/(parmDict[phfx+'Size;i']*parmDict[phfx+'Size;a']*costh)
     1124            Sgam *= np.sqrt((sinP*parmDict[phfx+'Size;a'])**2+(cosP*parmDict[phfx+'Size;i'])**2)
     1125        else:           #ellipsoidal crystallites
     1126            Sij =[parmDict[phfx+'Size:%d'%(i)] for i in range(6)]
     1127            H = np.array(refl[:3])
     1128            lenR = G2pwd.ellipseSize(H,Sij,GB)
     1129            Sgam = 1.8*wave/(np.pi*costh*lenR)
     1130        #microstrain               
     1131        if calcControls[phfx+'MustrainType'] == 'isotropic':
     1132            Mgam = 0.018*parmDict[phfx+'Mustrain;i']*tand(refl[5]/2.)/np.pi
     1133        elif calcControls[phfx+'MustrainType'] == 'uniaxial':
     1134            H = np.array(refl[:3])
     1135            P = np.array(calcControls[phfx+'MustrainAxis'])
     1136            cosP,sinP = G2lat.CosSinAngle(H,P,G)
     1137            Si = parmDict[phfx+'Mustrain;i']
     1138            Sa = parmDict[phfx+'Mustrain;a']
     1139            Mgam = 0.018*Si*Sa*tand(refl[5]/2.)/(np.pi*np.sqrt((Si*cosP)**2+(Sa*sinP)**2))
     1140        else:       #generalized - P.W. Stephens model
     1141            pwrs = calcControls[phfx+'MuPwrs']
     1142            sum = 0
     1143            for i,pwr in enumerate(pwrs):
     1144                sum += parmDict[phfx+'Mustrain:'+str(i)]*refl[0]**pwr[0]*refl[1]**pwr[1]*refl[2]**pwr[2]
     1145            Mgam = 0.018*refl[4]**2*tand(refl[5]/2.)*sum
     1146    elif 'T' in calcControls[hfx+'histType']:
     1147        #crystallite size
     1148        if calcControls[phfx+'SizeType'] == 'isotropic':
     1149            Sgam = 1.e-4*parmDict[hfx+'difC']*parmDict[phfx+'Size;i']
     1150        elif calcControls[phfx+'SizeType'] == 'uniaxial':
     1151            H = np.array(refl[:3])
     1152            P = np.array(calcControls[phfx+'SizeAxis'])
     1153            cosP,sinP = G2lat.CosSinAngle(H,P,G)
     1154            Sgam = 1.e-4*parmDict[hfx+'difC']*(parmDict[phfx+'Size;i']*parmDict[phfx+'Size;a'])
     1155            Sgam *= np.sqrt((sinP*parmDict[phfx+'Size;a'])**2+(cosP*parmDict[phfx+'Size;i'])**2)
     1156        else:           #ellipsoidal crystallites
     1157            Sij =[parmDict[phfx+'Size:%d'%(i)] for i in range(6)]
     1158            H = np.array(refl[:3])
     1159            lenR = G2pwd.ellipseSize(H,Sij,GB)
     1160            Sgam = 1.e-4*parmDict[hfx+'difC']*(refl[4]**2*lenR)
     1161        #microstrain               
     1162        if calcControls[phfx+'MustrainType'] == 'isotropic':
     1163            Mgam = 1.e-6*parmDict[hfx+'difC']*parmDict[phfx+'Mustrain;i']
     1164        elif calcControls[phfx+'MustrainType'] == 'uniaxial':
     1165            H = np.array(refl[:3])
     1166            P = np.array(calcControls[phfx+'MustrainAxis'])
     1167            cosP,sinP = G2lat.CosSinAngle(H,P,G)
     1168            Si = parmDict[phfx+'Mustrain;i']
     1169            Sa = parmDict[phfx+'Mustrain;a']
     1170            Mgam = 1.e-6*parmDict[hfx+'difC']*Si*Sa/np.sqrt((Si*cosP)**2+(Sa*sinP)**2)
     1171        else:       #generalized - P.W. Stephens model
     1172            pwrs = calcControls[phfx+'MuPwrs']
     1173            sum = 0
     1174            for i,pwr in enumerate(pwrs):
     1175                sum += parmDict[phfx+'Mustrain:'+str(i)]*refl[0]**pwr[0]*refl[1]**pwr[1]*refl[2]**pwr[2]
     1176            Mgam = 1.e-6*parmDict[hfx+'difC']*refl[4]**2*sum
     1177           
    11191178    gam = Sgam*parmDict[phfx+'Size;mx']+Mgam*parmDict[phfx+'Mustrain;mx']
    11201179    sig = (Sgam*(1.-parmDict[phfx+'Size;mx']))**2+(Mgam*(1.-parmDict[phfx+'Mustrain;mx']))**2
     
    11221181    return sig,gam
    11231182       
    1124 def GetSampleSigGamDerv(refl,wave,G,GB,phfx,calcControls,parmDict):
     1183def GetSampleSigGamDerv(refl,wave,G,GB,hfx,phfx,calcControls,parmDict):
    11251184    'Needs a doc string'
    11261185    gamDict = {}
    11271186    sigDict = {}
    1128     costh = cosd(refl[5]/2.)
    1129     tanth = tand(refl[5]/2.)
    1130     #crystallite size derivatives
    1131     if calcControls[phfx+'SizeType'] == 'isotropic':
    1132         Sgam = 1.8*wave/(np.pi*parmDict[phfx+'Size;i']*costh)
    1133         gamDict[phfx+'Size;i'] = -1.8*wave*parmDict[phfx+'Size;mx']/(np.pi*costh)
    1134         sigDict[phfx+'Size;i'] = -3.6*Sgam*wave*(1.-parmDict[phfx+'Size;mx'])**2/(np.pi*costh*ateln2)
    1135     elif calcControls[phfx+'SizeType'] == 'uniaxial':
    1136         H = np.array(refl[:3])
    1137         P = np.array(calcControls[phfx+'SizeAxis'])
    1138         cosP,sinP = G2lat.CosSinAngle(H,P,G)
    1139         Si = parmDict[phfx+'Size;i']
    1140         Sa = parmDict[phfx+'Size;a']
    1141         gami = (1.8*wave/np.pi)/(Si*Sa)
    1142         sqtrm = np.sqrt((sinP*Sa)**2+(cosP*Si)**2)
    1143         Sgam = gami*sqtrm
    1144         gam = Sgam/costh
    1145         dsi = (gami*Si*cosP**2/(sqtrm*costh)-gam/Si)
    1146         dsa = (gami*Sa*sinP**2/(sqtrm*costh)-gam/Sa)
    1147         gamDict[phfx+'Size;i'] = dsi*parmDict[phfx+'Size;mx']
    1148         gamDict[phfx+'Size;a'] = dsa*parmDict[phfx+'Size;mx']
    1149         sigDict[phfx+'Size;i'] = 2.*dsi*Sgam*(1.-parmDict[phfx+'Size;mx'])**2/ateln2
    1150         sigDict[phfx+'Size;a'] = 2.*dsa*Sgam*(1.-parmDict[phfx+'Size;mx'])**2/ateln2
    1151     else:           #ellipsoidal crystallites
    1152         const = 1.8*wave/(np.pi*costh)
    1153         Sij =[parmDict[phfx+'Size:%d'%(i)] for i in range(6)]
    1154         H = np.array(refl[:3])
    1155         lenR,dRdS = G2pwd.ellipseSizeDerv(H,Sij,GB)
    1156         Sgam = 1.8*wave/(np.pi*costh*lenR)
    1157         for i,item in enumerate([phfx+'Size:%d'%(j) for j in range(6)]):
    1158             gamDict[item] = -(const/lenR**2)*dRdS[i]*parmDict[phfx+'Size;mx']
    1159             sigDict[item] = -2.*Sgam*(const/lenR**2)*dRdS[i]*(1.-parmDict[phfx+'Size;mx'])**2/ateln2
    1160     gamDict[phfx+'Size;mx'] = Sgam
    1161     sigDict[phfx+'Size;mx'] = -2.*Sgam**2*(1.-parmDict[phfx+'Size;mx'])/ateln2
    1162            
    1163     #microstrain derivatives               
    1164     if calcControls[phfx+'MustrainType'] == 'isotropic':
    1165         Mgam = 0.018*parmDict[phfx+'Mustrain;i']*tand(refl[5]/2.)/np.pi
    1166         gamDict[phfx+'Mustrain;i'] =  0.018*tanth*parmDict[phfx+'Mustrain;mx']/np.pi
    1167         sigDict[phfx+'Mustrain;i'] =  0.036*Mgam*tanth*(1.-parmDict[phfx+'Mustrain;mx'])**2/(np.pi*ateln2)       
    1168     elif calcControls[phfx+'MustrainType'] == 'uniaxial':
    1169         H = np.array(refl[:3])
    1170         P = np.array(calcControls[phfx+'MustrainAxis'])
    1171         cosP,sinP = G2lat.CosSinAngle(H,P,G)
    1172         Si = parmDict[phfx+'Mustrain;i']
    1173         Sa = parmDict[phfx+'Mustrain;a']
    1174         gami = 0.018*Si*Sa*tanth/np.pi
    1175         sqtrm = np.sqrt((Si*cosP)**2+(Sa*sinP)**2)
    1176         Mgam = gami/sqtrm
    1177         dsi = -gami*Si*cosP**2/sqtrm**3
    1178         dsa = -gami*Sa*sinP**2/sqtrm**3
    1179         gamDict[phfx+'Mustrain;i'] = (Mgam/Si+dsi)*parmDict[phfx+'Mustrain;mx']
    1180         gamDict[phfx+'Mustrain;a'] = (Mgam/Sa+dsa)*parmDict[phfx+'Mustrain;mx']
    1181         sigDict[phfx+'Mustrain;i'] = 2*(Mgam/Si+dsi)*Mgam*(1.-parmDict[phfx+'Mustrain;mx'])**2/ateln2
    1182         sigDict[phfx+'Mustrain;a'] = 2*(Mgam/Sa+dsa)*Mgam*(1.-parmDict[phfx+'Mustrain;mx'])**2/ateln2       
    1183     else:       #generalized - P.W. Stephens model
    1184         pwrs = calcControls[phfx+'MuPwrs']
    1185         const = 0.018*refl[4]**2*tanth
    1186         sum = 0
    1187         for i,pwr in enumerate(pwrs):
    1188             term = refl[0]**pwr[0]*refl[1]**pwr[1]*refl[2]**pwr[2]
    1189             sum += parmDict[phfx+'Mustrain:'+str(i)]*term
    1190             gamDict[phfx+'Mustrain:'+str(i)] = const*term*parmDict[phfx+'Mustrain;mx']
    1191             sigDict[phfx+'Mustrain:'+str(i)] = \
    1192                 2.*const*term*(1.-parmDict[phfx+'Mustrain;mx'])**2/ateln2
    1193         Mgam = 0.018*refl[4]**2*tand(refl[5]/2.)*sum
    1194         for i in range(len(pwrs)):
    1195             sigDict[phfx+'Mustrain:'+str(i)] *= Mgam
    1196     gamDict[phfx+'Mustrain;mx'] = Mgam
    1197     sigDict[phfx+'Mustrain;mx'] = -2.*Mgam**2*(1.-parmDict[phfx+'Mustrain;mx'])/ateln2
     1187    if 'C' in calcControls[hfx+'histType']:
     1188        costh = cosd(refl[5]/2.)
     1189        tanth = tand(refl[5]/2.)
     1190        #crystallite size derivatives
     1191        if calcControls[phfx+'SizeType'] == 'isotropic':
     1192            Sgam = 1.8*wave/(np.pi*parmDict[phfx+'Size;i']*costh)
     1193            gamDict[phfx+'Size;i'] = -1.8*wave*parmDict[phfx+'Size;mx']/(np.pi*costh)
     1194            sigDict[phfx+'Size;i'] = -3.6*Sgam*wave*(1.-parmDict[phfx+'Size;mx'])**2/(np.pi*costh*ateln2)
     1195        elif calcControls[phfx+'SizeType'] == 'uniaxial':
     1196            H = np.array(refl[:3])
     1197            P = np.array(calcControls[phfx+'SizeAxis'])
     1198            cosP,sinP = G2lat.CosSinAngle(H,P,G)
     1199            Si = parmDict[phfx+'Size;i']
     1200            Sa = parmDict[phfx+'Size;a']
     1201            gami = (1.8*wave/np.pi)/(Si*Sa)
     1202            sqtrm = np.sqrt((sinP*Sa)**2+(cosP*Si)**2)
     1203            Sgam = gami*sqtrm
     1204            gam = Sgam/costh
     1205            dsi = (gami*Si*cosP**2/(sqtrm*costh)-gam/Si)
     1206            dsa = (gami*Sa*sinP**2/(sqtrm*costh)-gam/Sa)
     1207            gamDict[phfx+'Size;i'] = dsi*parmDict[phfx+'Size;mx']
     1208            gamDict[phfx+'Size;a'] = dsa*parmDict[phfx+'Size;mx']
     1209            sigDict[phfx+'Size;i'] = 2.*dsi*Sgam*(1.-parmDict[phfx+'Size;mx'])**2/ateln2
     1210            sigDict[phfx+'Size;a'] = 2.*dsa*Sgam*(1.-parmDict[phfx+'Size;mx'])**2/ateln2
     1211        else:           #ellipsoidal crystallites
     1212            const = 1.8*wave/(np.pi*costh)
     1213            Sij =[parmDict[phfx+'Size:%d'%(i)] for i in range(6)]
     1214            H = np.array(refl[:3])
     1215            lenR,dRdS = G2pwd.ellipseSizeDerv(H,Sij,GB)
     1216            Sgam = const/lenR
     1217            for i,item in enumerate([phfx+'Size:%d'%(j) for j in range(6)]):
     1218                gamDict[item] = -(const/lenR**2)*dRdS[i]*parmDict[phfx+'Size;mx']
     1219                sigDict[item] = -2.*Sgam*(const/lenR**2)*dRdS[i]*(1.-parmDict[phfx+'Size;mx'])**2/ateln2
     1220        gamDict[phfx+'Size;mx'] = Sgam
     1221        sigDict[phfx+'Size;mx'] = -2.*Sgam**2*(1.-parmDict[phfx+'Size;mx'])/ateln2
     1222               
     1223        #microstrain derivatives               
     1224        if calcControls[phfx+'MustrainType'] == 'isotropic':
     1225            Mgam = 0.018*parmDict[phfx+'Mustrain;i']*tand(refl[5]/2.)/np.pi
     1226            gamDict[phfx+'Mustrain;i'] =  0.018*tanth*parmDict[phfx+'Mustrain;mx']/np.pi
     1227            sigDict[phfx+'Mustrain;i'] =  0.036*Mgam*tanth*(1.-parmDict[phfx+'Mustrain;mx'])**2/(np.pi*ateln2)       
     1228        elif calcControls[phfx+'MustrainType'] == 'uniaxial':
     1229            H = np.array(refl[:3])
     1230            P = np.array(calcControls[phfx+'MustrainAxis'])
     1231            cosP,sinP = G2lat.CosSinAngle(H,P,G)
     1232            Si = parmDict[phfx+'Mustrain;i']
     1233            Sa = parmDict[phfx+'Mustrain;a']
     1234            gami = 0.018*Si*Sa*tanth/np.pi
     1235            sqtrm = np.sqrt((Si*cosP)**2+(Sa*sinP)**2)
     1236            Mgam = gami/sqtrm
     1237            dsi = -gami*Si*cosP**2/sqtrm**3
     1238            dsa = -gami*Sa*sinP**2/sqtrm**3
     1239            gamDict[phfx+'Mustrain;i'] = (Mgam/Si+dsi)*parmDict[phfx+'Mustrain;mx']
     1240            gamDict[phfx+'Mustrain;a'] = (Mgam/Sa+dsa)*parmDict[phfx+'Mustrain;mx']
     1241            sigDict[phfx+'Mustrain;i'] = 2*(Mgam/Si+dsi)*Mgam*(1.-parmDict[phfx+'Mustrain;mx'])**2/ateln2
     1242            sigDict[phfx+'Mustrain;a'] = 2*(Mgam/Sa+dsa)*Mgam*(1.-parmDict[phfx+'Mustrain;mx'])**2/ateln2       
     1243        else:       #generalized - P.W. Stephens model
     1244            pwrs = calcControls[phfx+'MuPwrs']
     1245            const = 0.018*refl[4]**2*tanth
     1246            sum = 0
     1247            for i,pwr in enumerate(pwrs):
     1248                term = refl[0]**pwr[0]*refl[1]**pwr[1]*refl[2]**pwr[2]
     1249                sum += parmDict[phfx+'Mustrain:'+str(i)]*term
     1250                gamDict[phfx+'Mustrain:'+str(i)] = const*term*parmDict[phfx+'Mustrain;mx']
     1251                sigDict[phfx+'Mustrain:'+str(i)] = \
     1252                    2.*const*term*(1.-parmDict[phfx+'Mustrain;mx'])**2/ateln2
     1253            Mgam = 0.018*refl[4]**2*tand(refl[5]/2.)*sum
     1254            for i in range(len(pwrs)):
     1255                sigDict[phfx+'Mustrain:'+str(i)] *= Mgam
     1256        gamDict[phfx+'Mustrain;mx'] = Mgam
     1257        sigDict[phfx+'Mustrain;mx'] = -2.*Mgam**2*(1.-parmDict[phfx+'Mustrain;mx'])/ateln2
     1258    else:   #'T'OF
     1259        if calcControls[phfx+'SizeType'] == 'isotropic':
     1260            Sgam = 1.e-4*parmDict[hfx+'difC']*parmDict[phfx+'Size;i']
     1261            gamDict[phfx+'Size;i'] = 1.e-4*parmDict[hfx+'difC']*parmDict[phfx+'Size;mx']
     1262            sigDict[phfx+'Size;i'] = 1.e-4*parmDict[hfx+'difC']*Sgam*(1.-parmDict[phfx+'Size;mx'])**2/ateln2
     1263        elif calcControls[phfx+'SizeType'] == 'uniaxial':
     1264            const = 1.e-4*parmDict[hfx+'difC']
     1265            H = np.array(refl[:3])
     1266            P = np.array(calcControls[phfx+'SizeAxis'])
     1267            cosP,sinP = G2lat.CosSinAngle(H,P,G)
     1268            Si = parmDict[phfx+'Size;i']
     1269            Sa = parmDict[phfx+'Size;a']
     1270            gami = const*(Si*Sa)
     1271            sqtrm = np.sqrt((sinP*Sa)**2+(cosP*Si)**2)
     1272            Sgam = gami*sqtrm
     1273            dsi = gami*Si*cosP**2/sqtrm-gam/Si
     1274            dsa = gami*Sa*sinP**2/sqtrm-gam/Sa
     1275            gamDict[phfx+'Size;i'] = const*parmDict[phfx+'Size;mx']*Sa
     1276            gamDict[phfx+'Size;a'] = const*parmDict[phfx+'Size;mx']*Si
     1277            sigDict[phfx+'Size;i'] = 2.*dsi*Sgam*(1.-parmDict[phfx+'Size;mx'])**2/ateln2
     1278            sigDict[phfx+'Size;a'] = 2.*dsa*Sgam*(1.-parmDict[phfx+'Size;mx'])**2/ateln2
     1279        else:           #ellipsoidal crystallites
     1280            const = 1.e-4*parmDict[hfx+'difC']
     1281            Sij =[parmDict[phfx+'Size:%d'%(i)] for i in range(6)]
     1282            H = np.array(refl[:3])
     1283            lenR,dRdS = G2pwd.ellipseSizeDerv(H,Sij,GB)
     1284            Sgam = const/lenR
     1285            for i,item in enumerate([phfx+'Size:%d'%(j) for j in range(6)]):
     1286                gamDict[item] = -(const/lenR**2)*dRdS[i]*parmDict[phfx+'Size;mx']
     1287                sigDict[item] = -2.*Sgam*(const/lenR**2)*dRdS[i]*(1.-parmDict[phfx+'Size;mx'])**2/ateln2
     1288        gamDict[phfx+'Size;mx'] = Sgam
     1289        sigDict[phfx+'Size;mx'] = -2.*Sgam**2*(1.-parmDict[phfx+'Size;mx'])/ateln2
     1290               
     1291        #microstrain derivatives               
     1292        if calcControls[phfx+'MustrainType'] == 'isotropic':
     1293            Mgam = 1.e-6*parmDict[hfx+'difC']*parmDict[phfx+'Mustrain;i']
     1294            gamDict[phfx+'Mustrain;i'] =  1.e-6*parmDict[hfx+'difC']*parmDict[phfx+'Mustrain;mx']
     1295            sigDict[phfx+'Mustrain;i'] =  2.e-6*parmDict[hfx+'difC']*Mgam*(1.-parmDict[phfx+'Mustrain;mx'])**2/ateln2       
     1296        elif calcControls[phfx+'MustrainType'] == 'uniaxial':
     1297            H = np.array(refl[:3])
     1298            P = np.array(calcControls[phfx+'MustrainAxis'])
     1299            cosP,sinP = G2lat.CosSinAngle(H,P,G)
     1300            Si = parmDict[phfx+'Mustrain;i']
     1301            Sa = parmDict[phfx+'Mustrain;a']
     1302            gami = 1.e-6*parmDict[hfx+'difC']*Si*Sa
     1303            sqtrm = np.sqrt((Si*cosP)**2+(Sa*sinP)**2)
     1304            Mgam = gami/sqtrm
     1305            dsi = -gami*Si*cosP**2/sqtrm**3
     1306            dsa = -gami*Sa*sinP**2/sqtrm**3
     1307            gamDict[phfx+'Mustrain;i'] = (Mgam/Si+dsi)*parmDict[phfx+'Mustrain;mx']
     1308            gamDict[phfx+'Mustrain;a'] = (Mgam/Sa+dsa)*parmDict[phfx+'Mustrain;mx']
     1309            sigDict[phfx+'Mustrain;i'] = 2*(Mgam/Si+dsi)*Mgam*(1.-parmDict[phfx+'Mustrain;mx'])**2/ateln2
     1310            sigDict[phfx+'Mustrain;a'] = 2*(Mgam/Sa+dsa)*Mgam*(1.-parmDict[phfx+'Mustrain;mx'])**2/ateln2       
     1311        else:       #generalized - P.W. Stephens model
     1312            pwrs = calcControls[phfx+'MuPwrs']
     1313            const = 1.e-6*parmDict[hfx+'difC']*refl[4]**2
     1314            sum = 0
     1315            for i,pwr in enumerate(pwrs):
     1316                term = refl[0]**pwr[0]*refl[1]**pwr[1]*refl[2]**pwr[2]
     1317                sum += parmDict[phfx+'Mustrain:'+str(i)]*term
     1318                gamDict[phfx+'Mustrain:'+str(i)] = const*term*parmDict[phfx+'Mustrain;mx']
     1319                sigDict[phfx+'Mustrain:'+str(i)] = \
     1320                    2.*const*term*(1.-parmDict[phfx+'Mustrain;mx'])**2/ateln2
     1321            Mgam = const*sum
     1322            for i in range(len(pwrs)):
     1323                sigDict[phfx+'Mustrain:'+str(i)] *= Mgam       
     1324        gamDict[phfx+'Mustrain;mx'] = Mgam
     1325        sigDict[phfx+'Mustrain;mx'] = -2.*Mgam**2*(1.-parmDict[phfx+'Mustrain;mx'])/ateln2
     1326       
    11981327    return sigDict,gamDict
    11991328       
     
    12051334
    12061335    refl[4] = d
    1207     pos = 2.0*asind(wave/(2.0*d))+parmDict[hfx+'Zero']
    1208     const = 9.e-2/(np.pi*parmDict[hfx+'Gonio. radius'])                  #shifts in microns
    1209     if 'Bragg' in calcControls[hfx+'instType']:
    1210         pos -= const*(4.*parmDict[hfx+'Shift']*cosd(pos/2.0)+ \
    1211             parmDict[hfx+'Transparency']*sind(pos)*100.0)            #trans(=1/mueff) in cm
    1212     else:               #Debye-Scherrer - simple but maybe not right
    1213         pos -= const*(parmDict[hfx+'DisplaceX']*cosd(pos)+parmDict[hfx+'DisplaceY']*sind(pos))
     1336    if 'C' in calcControls[hfx+'histType']:
     1337        pos = 2.0*asind(wave/(2.0*d))+parmDict[hfx+'Zero']
     1338        const = 9.e-2/(np.pi*parmDict[hfx+'Gonio. radius'])                  #shifts in microns
     1339        if 'Bragg' in calcControls[hfx+'instType']:
     1340            pos -= const*(4.*parmDict[hfx+'Shift']*cosd(pos/2.0)+ \
     1341                parmDict[hfx+'Transparency']*sind(pos)*100.0)            #trans(=1/mueff) in cm
     1342        else:               #Debye-Scherrer - simple but maybe not right
     1343            pos -= const*(parmDict[hfx+'DisplaceX']*cosd(pos)+parmDict[hfx+'DisplaceY']*sind(pos))
     1344    elif 'T' in calcControls[hfx+'histType']:
     1345        pos = parmDict[hfx+'difC']*d+parmDict[hfx+'difA']*d**2+parmDict[hfx+'difB']*d**3+parmDict[hfx+'Zero']
     1346        #do I need sample position effects - maybe?
    12141347    return pos
    12151348
     
    12201353    dstsq = G2lat.calc_rDsq(np.array([h,k,l]),A)
    12211354    dst = np.sqrt(dstsq)
    1222     pos = refl[5]-parmDict[hfx+'Zero']
    1223     const = dpr/np.sqrt(1.0-wave**2*dstsq/4.0)
    1224     dpdw = const*dst
    1225     dpdA = np.array([h**2,k**2,l**2,h*k,h*l,k*l])
    1226     dpdA *= const*wave/(2.0*dst)
    1227     dpdZ = 1.0
    1228     const = 9.e-2/(np.pi*parmDict[hfx+'Gonio. radius'])                  #shifts in microns
    1229     if 'Bragg' in calcControls[hfx+'instType']:
    1230         dpdSh = -4.*const*cosd(pos/2.0)
    1231         dpdTr = -const*sind(pos)*100.0
    1232         return dpdA,dpdw,dpdZ,dpdSh,dpdTr,0.,0.
    1233     else:               #Debye-Scherrer - simple but maybe not right
    1234         dpdXd = -const*cosd(pos)
    1235         dpdYd = -const*sind(pos)
    1236         return dpdA,dpdw,dpdZ,0.,0.,dpdXd,dpdYd
     1355    dsp = 1./dst
     1356    if 'C' in calcControls[hfx+'histType']:
     1357        pos = refl[5]-parmDict[hfx+'Zero']
     1358        const = dpr/np.sqrt(1.0-wave**2*dstsq/4.0)
     1359        dpdw = const*dst
     1360        dpdA = np.array([h**2,k**2,l**2,h*k,h*l,k*l])
     1361        dpdA *= const*wave/(2.0*dst)
     1362        dpdZ = 1.0
     1363        const = 9.e-2/(np.pi*parmDict[hfx+'Gonio. radius'])                  #shifts in microns
     1364        if 'Bragg' in calcControls[hfx+'instType']:
     1365            dpdSh = -4.*const*cosd(pos/2.0)
     1366            dpdTr = -const*sind(pos)*100.0
     1367            return dpdA,dpdw,dpdZ,dpdSh,dpdTr,0.,0.
     1368        else:               #Debye-Scherrer - simple but maybe not right
     1369            dpdXd = -const*cosd(pos)
     1370            dpdYd = -const*sind(pos)
     1371            return dpdA,dpdw,dpdZ,0.,0.,dpdXd,dpdYd
     1372    elif 'T' in calcControls[hfx+'histType']:
     1373        dpdA = np.array([h**2,k**2,l**2,h*k,h*l,k*l])
     1374        dpdZ = 1.0
     1375        dpdDC = dsp
     1376        dpdDA = dsp**2
     1377        dpdDB = dsp**3
     1378        return dpdA,dpdZ,dpdDC,dpdDA,dpdDB
    12371379           
    12381380def GetHStrainShift(refl,SGData,phfx,parmDict):
     
    13071449            hfx = ':%d:'%(hId)
    13081450            Limits = calcControls[hfx+'Limits']
    1309             shl = max(parmDict[hfx+'SH/L'],0.0005)
    1310             Ka2 = False
    1311             kRatio = 0.0
    1312             if hfx+'Lam1' in parmDict.keys():
    1313                 Ka2 = True
    1314                 lamRatio = 360*(parmDict[hfx+'Lam2']-parmDict[hfx+'Lam1'])/(np.pi*parmDict[hfx+'Lam1'])
    1315                 kRatio = parmDict[hfx+'I(L2)/I(L1)']
     1451            if 'C' in calcControls[hfx+'histType']:
     1452                shl = max(parmDict[hfx+'SH/L'],0.0005)
     1453                Ka2 = False
     1454                kRatio = 0.0
     1455                if hfx+'Lam1' in parmDict.keys():
     1456                    Ka2 = True
     1457                    lamRatio = 360*(parmDict[hfx+'Lam2']-parmDict[hfx+'Lam1'])/(np.pi*parmDict[hfx+'Lam1'])
     1458                    kRatio = parmDict[hfx+'I(L2)/I(L1)']
    13161459            x,y,w,yc,yb,yd = Histogram['Data']
    13171460            xB = np.searchsorted(x,Limits[0])
     
    13381481                        iFin = max(xB,min(np.searchsorted(x,refl[5]+fmax),xF))
    13391482                        iFin2 = iFin
    1340                         yp[iBeg:iFin] = refl[11]*refl[9]*G2pwd.getFCJVoigt3(refl[5],refl[6],refl[7],shl,x[iBeg:iFin])    #>90% of time spent here
     1483                        yp[iBeg:iFin] = refl[11]*refl[9]*G2pwd.getFCJVoigt3(refl[5],refl[6],refl[7],shl,ma.getdata(x[iBeg:iFin]))    #>90% of time spent here
    13411484                        if Ka2:
    13421485                            pos2 = refl[5]+lamRatio*tand(refl[5]/2.0)       # + 360/pi * Dlam/lam * tan(th)
     
    13481491                            elif not iBeg2-iFin2:     #peak above high limit - done
    13491492                                break
    1350                             yp[iBeg2:iFin2] += refl[11]*refl[9]*kRatio*G2pwd.getFCJVoigt3(pos2,refl[6],refl[7],shl,x[iBeg2:iFin2])        #and here
     1493                            yp[iBeg2:iFin2] += refl[11]*refl[9]*kRatio*G2pwd.getFCJVoigt3(pos2,refl[6],refl[7],shl,ma.getdata(x[iBeg2:iFin2]))        #and here
    13511494                        refl[8] = np.sum(np.where(ratio[iBeg:iFin2]>0.,yp[iBeg:iFin2]*ratio[iBeg:iFin2]/(refl[11]*(1.+kRatio)),0.0))
    13521495                    elif 'T' in calcControls[hfx+'histType']:
    1353                         print 'TOF Undefined at present'
    1354                         raise Exception    #no TOF yet
     1496                        yp = np.zeros_like(yb)
     1497                        Wd,fmin,fmax = G2pwd.getWidthsTOF(refl[5],refl[12],refl[13],refl[6],refl[7])
     1498                        iBeg = max(xB,np.searchsorted(x,refl[5]-fmin))
     1499                        iFin = max(xB,min(np.searchsorted(x,refl[5]+fmax),xF))
     1500                        yp[iBeg:iFin] = refl[11]*refl[9]*G2pwd.getEpsVoigt(refl[5],refl[12],refl[13],refl[6],refl[7],ma.getdata(x[iBeg:iFin]))  #>90% of time spent here
     1501                        refl[8] = np.sum(np.where(ratio[iBeg:iFin]>0.,yp[iBeg:iFin]*ratio[iBeg:iFin]/refl[11],0.0))
    13551502                    Fo = np.sqrt(np.abs(refl[8]))
    13561503                    Fc = np.sqrt(np.abs(refl[9]))
     
    13701517    'Needs a doc string'
    13711518   
    1372     def GetReflSigGam(refl,wave,G,GB,hfx,phfx,calcControls,parmDict):
     1519    def GetReflSigGamCW(refl,wave,G,GB,phfx,calcControls,parmDict):
    13731520        U = parmDict[hfx+'U']
    13741521        V = parmDict[hfx+'V']
     
    13771524        Y = parmDict[hfx+'Y']
    13781525        tanPos = tand(refl[5]/2.0)
    1379         Ssig,Sgam = GetSampleSigGam(refl,wave,G,GB,phfx,calcControls,parmDict)
     1526        Ssig,Sgam = GetSampleSigGam(refl,wave,G,GB,hfx,phfx,calcControls,parmDict)
    13801527        sig = U*tanPos**2+V*tanPos+W+Ssig     #save peak sigma
    13811528        sig = max(0.001,sig)
     
    13831530        gam = max(0.001,gam)
    13841531        return sig,gam
     1532               
     1533    def GetReflSigGamTOF(refl,G,GB,phfx,calcControls,parmDict):
     1534        sig = parmDict[hfx+'sig-0']+parmDict[hfx+'sig-1']*refl[4]**2+parmDict[hfx+'sig-q']*refl[4]
     1535        gam = parmDict[hfx+'X']*refl[4]+parmDict[hfx+'Y']*refl[4]**2
     1536        Ssig,Sgam = GetSampleSigGam(refl,0.0,G,GB,hfx,phfx,calcControls,parmDict)
     1537        sig += Ssig     #save peak sigma
     1538        sig = max(0.001,sig)
     1539        gam += Sgam     #save peak gamma
     1540        gam = max(0.001,gam)
     1541        return sig,gam
     1542       
     1543    def GetReflAlpBet(refl,hfx,parmDict):
     1544        alp = parmDict[hfx+'alpha']/refl[4]
     1545        bet = parmDict[hfx+'beta-0']+parmDict[hfx+'beta-1']/refl[4]**4+parmDict[hfx+'beta-q']/refl[4]
     1546        return alp,bet
    13851547               
    13861548    hId = Histogram['hId']
     
    14001562        else:
    14011563            wave = parmDict[hfx+'Lam']
    1402     else:
    1403         print 'TOF Undefined at present - might be nothing need be done here'
    14041564    for phase in Histogram['Reflection Lists']:
    14051565        refDict = Histogram['Reflection Lists'][phase]
     
    14271587                Lorenz = 1./(2.*sind(refl[5]/2.)**2*cosd(refl[5]/2.))           #Lorentz correction
    14281588                refl[5] += GetHStrainShift(refl,SGData,phfx,parmDict)               #apply hydrostatic strain shift
    1429                 refl[6:8] = GetReflSigGam(refl,wave,G,GB,hfx,phfx,calcControls,parmDict)    #peak sig & gam
    1430                 GetIntensityCorr(refl,Uniq,G,g,pfx,phfx,hfx,SGData,calcControls,parmDict)    #puts corrections in refl[11]
    1431                 refl[11] *= Vst*Lorenz
     1589                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
     1591                 
    14321592                if Phase['General'].get('doPawley'):
    14331593                    try:
     
    14561616                    yc[iBeg:iFin] += refl[11]*refl[9]*kRatio*G2pwd.getFCJVoigt3(pos2,refl[6],refl[7],shl,ma.getdata(x[iBeg:iFin]))        #and here
    14571617            elif 'T' in calcControls[hfx+'histType']:
    1458                 print 'TOF Undefined at present'
    1459                 raise Exception    #no TOF yet
     1618                h,k,l = refl[:3]
     1619                Uniq = np.inner(refl[:3],SGMT)
     1620                refl[5] = GetReflPos(refl,0.0,G,hfx,calcControls,parmDict)         #corrected reflection position
     1621                Lorenz = sind(parmDict[hfx+'2-theta']/2)*refl[4]**4                                                #TOF Lorentz correction
     1622                refl[5] += GetHStrainShift(refl,SGData,phfx,parmDict)               #apply hydrostatic strain shift
     1623                refl[6:8] = GetReflSigGamTOF(refl,G,GB,phfx,calcControls,parmDict)    #peak sig & gam
     1624                refl[12:14] = GetReflAlpBet(refl,hfx,parmDict)
     1625                refl[11] = GetIntensityCorr(refl,Uniq,G,g,pfx,phfx,hfx,SGData,calcControls,parmDict)*Vst*Lorenz
     1626                if Phase['General'].get('doPawley'):
     1627                    try:
     1628                        pInd =pfx+'PWLref:%d'%(pawleyLookup[pfx+'%d,%d,%d'%(h,k,l)])
     1629                        refl[9] = parmDict[pInd]
     1630                    except KeyError:
     1631#                        print ' ***Error %d,%d,%d missing from Pawley reflection list ***'%(h,k,l)
     1632                        continue
     1633                Wd,fmin,fmax = G2pwd.getWidthsTOF(refl[5],refl[12],refl[13],refl[6],refl[7])
     1634                iBeg = np.searchsorted(x,refl[5]-fmin)
     1635                iFin = np.searchsorted(x,refl[5]+fmax)
     1636                if not iBeg+iFin:       #peak below low limit - skip peak
     1637                    continue
     1638                elif not iBeg-iFin:     #peak above high limit - done
     1639                    break
     1640                yc[iBeg:iFin] += refl[11]*refl[9]*G2pwd.getEpsVoigt(refl[5],refl[12],refl[13],refl[6],refl[7],ma.getdata(x[iBeg:iFin]))
    14601641#        print 'profile calc time: %.3fs'%(time.time()-time0)
    14611642    return yc,yb
     
    15181699    cw = np.diff(x)
    15191700    cw = np.append(cw,cw[-1])
     1701    Ka2 = False #also for TOF!
    15201702    if 'C' in calcControls[hfx+'histType']:   
    15211703        shl = max(parmDict[hfx+'SH/L'],0.002)
    1522         Ka2 = False
    15231704        if hfx+'Lam1' in parmDict.keys():
    15241705            wave = parmDict[hfx+'Lam1']
     
    15281709        else:
    15291710            wave = parmDict[hfx+'Lam']
    1530     else:
    1531         print 'TOF Undefined at present'
    1532         raise ValueError
    15331711    for phase in Histogram['Reflection Lists']:
    15341712        refDict = Histogram['Reflection Lists'][phase]
     
    15491727        time0 = time.time()
    15501728        for iref,refl in enumerate(refDict['RefList']):
     1729            h,k,l = refl[:3]
     1730            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)
    15511732            if 'C' in calcControls[hfx+'histType']:        #CW powder
    1552                 h,k,l = refl[:3]
    1553                 Uniq = np.inner(refl[:3],SGMT)
    1554                 dIdsh,dIdsp,dIdpola,dIdPO,dFdODF,dFdSA,dFdAb = GetIntensityDerv(refl,Uniq,G,g,pfx,phfx,hfx,SGData,calcControls,parmDict)
    15551733                Wd,fmin,fmax = G2pwd.getWidthsCW(refl[5],refl[6],refl[7],shl)
    1556                 iBeg = np.searchsorted(x,refl[5]-fmin)
    1557                 iFin = np.searchsorted(x,refl[5]+fmax)
    1558                 if not iBeg+iFin:       #peak below low limit - skip peak
    1559                     continue
    1560                 elif not iBeg-iFin:     #peak above high limit - done
    1561                     break
    1562                 pos = refl[5]
     1734            else: #'T'OF
     1735                Wd,fmin,fmax = G2pwd.getWidthsTOF(refl[5],refl[12],refl[13],refl[6],refl[7])
     1736            iBeg = np.searchsorted(x,refl[5]-fmin)
     1737            iFin = np.searchsorted(x,refl[5]+fmax)
     1738            if not iBeg+iFin:       #peak below low limit - skip peak
     1739                continue
     1740            elif not iBeg-iFin:     #peak above high limit - done
     1741                break
     1742            pos = refl[5]
     1743            if 'C' in calcControls[hfx+'histType']:
    15631744                tanth = tand(pos/2.0)
    15641745                costh = cosd(pos/2.0)
     
    15811762                        dMdpk2[5] = 100.*cw[iBeg2:iFin2]*refl[11]*dMdipk2[0]
    15821763                        dervDict2 = {'int':dMdpk2[0],'pos':dMdpk2[1],'sig':dMdpk2[2],'gam':dMdpk2[3],'shl':dMdpk2[4],'L1/L2':dMdpk2[5]*refl[9]}
    1583                 if Phase['General'].get('doPawley'):
    1584                     dMdpw = np.zeros(len(x))
    1585                     try:
    1586                         pIdx = pfx+'PWLref:'+str(pawleyLookup[pfx+'%d,%d,%d'%(h,k,l)])
    1587                         idx = varylist.index(pIdx)
    1588                         dMdpw[iBeg:iFin] = dervDict['int']/refl[9]
    1589                         if Ka2:
    1590                             dMdpw[iBeg2:iFin2] += dervDict2['int']/refl[9]
    1591                         dMdv[idx] = dMdpw
    1592                     except: # ValueError:
    1593                         pass
     1764            else:   #'T'OF
     1765                lenBF = iFin-iBeg
     1766                dMdpk = np.zeros(shape=(6,lenBF))
     1767                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]
     1770                dervDict = {'int':dMdpk[0],'pos':dMdpk[1],'alp':dMdpk[2],'bet':dMdpk[3],'sig':dMdpk[4],'gam':dMdpk[5]}           
     1771            if Phase['General'].get('doPawley'):
     1772                dMdpw = np.zeros(len(x))
     1773                try:
     1774                    pIdx = pfx+'PWLref:'+str(pawleyLookup[pfx+'%d,%d,%d'%(h,k,l)])
     1775                    idx = varylist.index(pIdx)
     1776                    dMdpw[iBeg:iFin] = dervDict['int']/refl[9]
     1777                    if Ka2: #not for TOF either
     1778                        dMdpw[iBeg2:iFin2] += dervDict2['int']/refl[9]
     1779                    dMdv[idx] = dMdpw
     1780                except: # ValueError:
     1781                    pass
     1782            if 'C' in calcControls[hfx+'histType']:
    15941783                dpdA,dpdw,dpdZ,dpdSh,dpdTr,dpdX,dpdY = GetReflPosDerv(refl,wave,A,hfx,calcControls,parmDict)
    15951784                names = {hfx+'Scale':[dIdsh,'int'],hfx+'Polariz.':[dIdpola,'int'],phfx+'Scale':[dIdsp,'int'],
     
    16041793                else:
    16051794                    names.update({hfx+'Absorption':[dFdAb,'int'],})
    1606                 for name in names:
    1607                     item = names[name]
    1608                     if name in varylist:
    1609                         dMdv[varylist.index(name)][iBeg:iFin] += item[0]*dervDict[item[1]]
     1795            else:   #'T'OF
     1796                dpdA,dpdZ,dpdDC,dpdDA,dpdDB = GetReflPosDerv(refl,0.0,A,hfx,calcControls,parmDict)
     1797                names = {hfx+'Scale':[dIdsh,'int'],phfx+'Scale':[dIdsp,'int'],
     1798                    hfx+'difC':[dpdDC,'pos'],hfx+'difA':[dpdDA,'pos'],hfx+'difB':[dpdDB,'pos'],
     1799                    hfx+'Zero':[dpdZ,'pos'],hfx+'X':[refl[4],'gam'],hfx+'Y':[refl[4]**2,'gam'],
     1800                    hfx+'alpha':[1./refl[4],'alp'],hfx+'beta-0':[1.0,'bet'],hfx+'beta-1':[1./refl[4]**4,'bet'],
     1801                    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'],}
     1803            for name in names:
     1804                item = names[name]
     1805                if name in varylist:
     1806                    dMdv[varylist.index(name)][iBeg:iFin] += item[0]*dervDict[item[1]]
     1807                    if Ka2:
     1808                        dMdv[varylist.index(name)][iBeg2:iFin2] += item[0]*dervDict2[item[1]]
     1809                elif name in dependentVars:
     1810                    depDerivDict[name][iBeg:iFin] += item[0]*dervDict[item[1]]
     1811                    if Ka2:
     1812                        depDerivDict[name][iBeg2:iFin2] += item[0]*dervDict2[item[1]]
     1813            for iPO in dIdPO:
     1814                if iPO in varylist:
     1815                    dMdv[varylist.index(iPO)][iBeg:iFin] += dIdPO[iPO]*dervDict['int']
     1816                    if Ka2:
     1817                        dMdv[varylist.index(iPO)][iBeg2:iFin2] += dIdPO[iPO]*dervDict2['int']
     1818                elif iPO in dependentVars:
     1819                    depDerivDict[iPO][iBeg:iFin] += dIdPO[iPO]*dervDict['int']
     1820                    if Ka2:
     1821                        depDerivDict[iPO][iBeg2:iFin2] += dIdPO[iPO]*dervDict2['int']
     1822            for i,name in enumerate(['omega','chi','phi']):
     1823                aname = pfx+'SH '+name
     1824                if aname in varylist:
     1825                    dMdv[varylist.index(aname)][iBeg:iFin] += dFdSA[i]*dervDict['int']
     1826                    if Ka2:
     1827                        dMdv[varylist.index(aname)][iBeg2:iFin2] += dFdSA[i]*dervDict2['int']
     1828                elif aname in dependentVars:
     1829                    depDerivDict[aname][iBeg:iFin] += dFdSA[i]*dervDict['int']
     1830                    if Ka2:
     1831                        depDerivDict[aname][iBeg2:iFin2] += dFdSA[i]*dervDict2['int']
     1832            for iSH in dFdODF:
     1833                if iSH in varylist:
     1834                    dMdv[varylist.index(iSH)][iBeg:iFin] += dFdODF[iSH]*dervDict['int']
     1835                    if Ka2:
     1836                        dMdv[varylist.index(iSH)][iBeg2:iFin2] += dFdODF[iSH]*dervDict2['int']
     1837                elif iSH in dependentVars:
     1838                    depDerivDict[iSH][iBeg:iFin] += dFdODF[iSH]*dervDict['int']
     1839                    if Ka2:
     1840                        depDerivDict[iSH][iBeg2:iFin2] += dFdODF[iSH]*dervDict2['int']
     1841            cellDervNames = cellVaryDerv(pfx,SGData,dpdA)
     1842            for name,dpdA in cellDervNames:
     1843                if name in varylist:
     1844                    dMdv[varylist.index(name)][iBeg:iFin] += dpdA*dervDict['pos']
     1845                    if Ka2:
     1846                        dMdv[varylist.index(name)][iBeg2:iFin2] += dpdA*dervDict2['pos']
     1847                elif name in dependentVars:
     1848                    depDerivDict[name][iBeg:iFin] += dpdA*dervDict['pos']
     1849                    if Ka2:
     1850                        depDerivDict[name][iBeg2:iFin2] += dpdA*dervDict2['pos']
     1851            dDijDict = GetHStrainShiftDerv(refl,SGData,phfx)
     1852            for name in dDijDict:
     1853                if name in varylist:
     1854                    dMdv[varylist.index(name)][iBeg:iFin] += dDijDict[name]*dervDict['pos']
     1855                    if Ka2:
     1856                        dMdv[varylist.index(name)][iBeg2:iFin2] += dDijDict[name]*dervDict2['pos']
     1857                elif name in dependentVars:
     1858                    depDerivDict[name][iBeg:iFin] += dDijDict[name]*dervDict['pos']
     1859                    if Ka2:
     1860                        depDerivDict[name][iBeg2:iFin2] += dDijDict[name]*dervDict2['pos']
     1861            if 'C' in calcControls[hfx+'histType']:
     1862                sigDict,gamDict = GetSampleSigGamDerv(refl,wave,G,GB,hfx,phfx,calcControls,parmDict)
     1863            else:   #'T'OF
     1864                sigDict,gamDict = GetSampleSigGamDerv(refl,0.0,G,GB,hfx,phfx,calcControls,parmDict)
     1865            for name in gamDict:
     1866                if name in varylist:
     1867                    dMdv[varylist.index(name)][iBeg:iFin] += gamDict[name]*dervDict['gam']
     1868                    if Ka2:
     1869                        dMdv[varylist.index(name)][iBeg2:iFin2] += gamDict[name]*dervDict2['gam']
     1870                elif name in dependentVars:
     1871                    depDerivDict[name][iBeg:iFin] += gamDict[name]*dervDict['gam']
     1872                    if Ka2:
     1873                        depDerivDict[name][iBeg2:iFin2] += gamDict[name]*dervDict2['gam']
     1874            for name in sigDict:
     1875                if name in varylist:
     1876                    dMdv[varylist.index(name)][iBeg:iFin] += sigDict[name]*dervDict['sig']
     1877                    if Ka2:
     1878                        dMdv[varylist.index(name)][iBeg2:iFin2] += sigDict[name]*dervDict2['sig']
     1879                elif name in dependentVars:
     1880                    depDerivDict[name][iBeg:iFin] += sigDict[name]*dervDict['sig']
     1881                    if Ka2:
     1882                        depDerivDict[name][iBeg2:iFin2] += sigDict[name]*dervDict2['sig']
     1883            for name in ['BabA','BabU']:
     1884                if refl[9]:
     1885                    if phfx+name in varylist:
     1886                        dMdv[varylist.index(phfx+name)][iBeg:iFin] += dFdvDict[pfx+name][iref]*dervDict['int']/refl[9]
    16101887                        if Ka2:
    1611                             dMdv[varylist.index(name)][iBeg2:iFin2] += item[0]*dervDict2[item[1]]
    1612                     elif name in dependentVars:
    1613                         depDerivDict[name][iBeg:iFin] += item[0]*dervDict[item[1]]
     1888                            dMdv[varylist.index(phfx+name)][iBeg2:iFin2] += dFdvDict[pfx+name][iref]*dervDict2['int']/refl[9]
     1889                    elif phfx+name in dependentVars:                   
     1890                        depDerivDict[phfx+name][iBeg:iFin] += dFdvDict[pfx+name][iref]*dervDict['int']/refl[9]
    16141891                        if Ka2:
    1615                             depDerivDict[name][iBeg2:iFin2] += item[0]*dervDict2[item[1]]
    1616                 for iPO in dIdPO:
    1617                     if iPO in varylist:
    1618                         dMdv[varylist.index(iPO)][iBeg:iFin] += dIdPO[iPO]*dervDict['int']
    1619                         if Ka2:
    1620                             dMdv[varylist.index(iPO)][iBeg2:iFin2] += dIdPO[iPO]*dervDict2['int']
    1621                     elif iPO in dependentVars:
    1622                         depDerivDict[iPO][iBeg:iFin] += dIdPO[iPO]*dervDict['int']
    1623                         if Ka2:
    1624                             depDerivDict[iPO][iBeg2:iFin2] += dIdPO[iPO]*dervDict2['int']
    1625                 for i,name in enumerate(['omega','chi','phi']):
    1626                     aname = pfx+'SH '+name
    1627                     if aname in varylist:
    1628                         dMdv[varylist.index(aname)][iBeg:iFin] += dFdSA[i]*dervDict['int']
    1629                         if Ka2:
    1630                             dMdv[varylist.index(aname)][iBeg2:iFin2] += dFdSA[i]*dervDict2['int']
    1631                     elif aname in dependentVars:
    1632                         depDerivDict[aname][iBeg:iFin] += dFdSA[i]*dervDict['int']
    1633                         if Ka2:
    1634                             depDerivDict[aname][iBeg2:iFin2] += dFdSA[i]*dervDict2['int']
    1635                 for iSH in dFdODF:
    1636                     if iSH in varylist:
    1637                         dMdv[varylist.index(iSH)][iBeg:iFin] += dFdODF[iSH]*dervDict['int']
    1638                         if Ka2:
    1639                             dMdv[varylist.index(iSH)][iBeg2:iFin2] += dFdODF[iSH]*dervDict2['int']
    1640                     elif iSH in dependentVars:
    1641                         depDerivDict[iSH][iBeg:iFin] += dFdODF[iSH]*dervDict['int']
    1642                         if Ka2:
    1643                             depDerivDict[iSH][iBeg2:iFin2] += dFdODF[iSH]*dervDict2['int']
    1644                 cellDervNames = cellVaryDerv(pfx,SGData,dpdA)
    1645                 for name,dpdA in cellDervNames:
    1646                     if name in varylist:
    1647                         dMdv[varylist.index(name)][iBeg:iFin] += dpdA*dervDict['pos']
    1648                         if Ka2:
    1649                             dMdv[varylist.index(name)][iBeg2:iFin2] += dpdA*dervDict2['pos']
    1650                     elif name in dependentVars:
    1651                         depDerivDict[name][iBeg:iFin] += dpdA*dervDict['pos']
    1652                         if Ka2:
    1653                             depDerivDict[name][iBeg2:iFin2] += dpdA*dervDict2['pos']
    1654                 dDijDict = GetHStrainShiftDerv(refl,SGData,phfx)
    1655                 for name in dDijDict:
    1656                     if name in varylist:
    1657                         dMdv[varylist.index(name)][iBeg:iFin] += dDijDict[name]*dervDict['pos']
    1658                         if Ka2:
    1659                             dMdv[varylist.index(name)][iBeg2:iFin2] += dDijDict[name]*dervDict2['pos']
    1660                     elif name in dependentVars:
    1661                         depDerivDict[name][iBeg:iFin] += dDijDict[name]*dervDict['pos']
    1662                         if Ka2:
    1663                             depDerivDict[name][iBeg2:iFin2] += dDijDict[name]*dervDict2['pos']
    1664                 sigDict,gamDict = GetSampleSigGamDerv(refl,wave,G,GB,phfx,calcControls,parmDict)
    1665                 for name in gamDict:
    1666                     if name in varylist:
    1667                         dMdv[varylist.index(name)][iBeg:iFin] += gamDict[name]*dervDict['gam']
    1668                         if Ka2:
    1669                             dMdv[varylist.index(name)][iBeg2:iFin2] += gamDict[name]*dervDict2['gam']
    1670                     elif name in dependentVars:
    1671                         depDerivDict[name][iBeg:iFin] += gamDict[name]*dervDict['gam']
    1672                         if Ka2:
    1673                             depDerivDict[name][iBeg2:iFin2] += gamDict[name]*dervDict2['gam']
    1674                 for name in sigDict:
    1675                     if name in varylist:
    1676                         dMdv[varylist.index(name)][iBeg:iFin] += sigDict[name]*dervDict['sig']
    1677                         if Ka2:
    1678                             dMdv[varylist.index(name)][iBeg2:iFin2] += sigDict[name]*dervDict2['sig']
    1679                     elif name in dependentVars:
    1680                         depDerivDict[name][iBeg:iFin] += sigDict[name]*dervDict['sig']
    1681                         if Ka2:
    1682                             depDerivDict[name][iBeg2:iFin2] += sigDict[name]*dervDict2['sig']
    1683                 for name in ['BabA','BabU']:
    1684                     if refl[9]:
    1685                         if phfx+name in varylist:
    1686                             dMdv[varylist.index(phfx+name)][iBeg:iFin] += dFdvDict[pfx+name][iref]*dervDict['int']/refl[9]
    1687                             if Ka2:
    1688                                 dMdv[varylist.index(phfx+name)][iBeg2:iFin2] += dFdvDict[pfx+name][iref]*dervDict2['int']/refl[9]
    1689                         elif phfx+name in dependentVars:                   
    1690                             depDerivDict[phfx+name][iBeg:iFin] += dFdvDict[pfx+name][iref]*dervDict['int']/refl[9]
    1691                             if Ka2:
    1692                                 depDerivDict[phfx+name][iBeg2:iFin2] += dFdvDict[pfx+name][iref]*dervDict2['int']/refl[9]                 
    1693             elif 'T' in calcControls[hfx+'histType']:
    1694                 print 'TOF Undefined at present'
    1695                 raise Exception    #no TOF yet
     1892                            depDerivDict[phfx+name][iBeg2:iFin2] += dFdvDict[pfx+name][iref]*dervDict2['int']/refl[9]                 
    16961893            if not Phase['General'].get('doPawley'):
    16971894                #do atom derivatives -  for RB,F,X & U so far             
     
    17161913                        if Ka2:
    17171914                            depDerivDict[name][iBeg2:iFin2] += dFdvDict[name][iref]*corr2
    1718 #        print 'profile derv time: %.3fs'%(time.time()-time0)
     1915    #        print 'profile derv time: %.3fs'%(time.time()-time0)
    17191916    # now process derivatives in constraints
    17201917    G2mv.Dict2Deriv(varylist,depDerivDict,dMdv)
     
    17221919   
    17231920def dervHKLF(Histogram,Phase,calcControls,varylist,parmDict,rigidbodyDict):
    1724     '''Loop over reflections ina HKLF histogram and compute derivatives of the fitting
     1921    '''Loop over reflections in a HKLF histogram and compute derivatives of the fitting
    17251922    model (M) with respect to all parameters.  Independent and dependant dM/dp arrays
    17261923    are returned to either dervRefine or HessRefine.
     
    17541951                    for j,var in enumerate(varylist):
    17551952                        if var in dFdvDict:
    1756                             dMdvh[j][iref] = w*dFdvDict[var][iref]*dervCor*parmDict[phfx+'Scale']
     1953                            dMdvh[j][iref] = w*dFdvDict[var][iref]*parmDict[phfx+'Scale']*dervCor*ref[11]
    17571954                    for var in dependentVars:
    17581955                        if var in dFdvDict:
    1759                             depDerivDict[var][iref] = w*dFdvDict[var][iref]*dervCor*parmDict[phfx+'Scale']
     1956                            depDerivDict[var][iref] = w*dFdvDict[var][iref]*parmDict[phfx+'Scale']*dervCor*ref[11]
    17601957                    if phfx+'Scale' in varylist:
    17611958                        dMdvh[varylist.index(phfx+'Scale')][iref] = w*ref[9]*dervCor
     
    17831980                    for j,var in enumerate(varylist):
    17841981                        if var in dFdvDict:
    1785                             dMdvh[j][iref] = w*dFdvDict[var][iref]*dervCor*parmDict[phfx+'Scale']
     1982                            dMdvh[j][iref] = w*dFdvDict[var][iref]*parmDict[phfx+'Scale']*dervCor*ref[11]
    17861983                    for var in dependentVars:
    17871984                        if var in dFdvDict:
    1788                             depDerivDict[var][iref] = w*dFdvDict[var][iref]*dervCor*parmDict[phfx+'Scale']
     1985                            depDerivDict[var][iref] = w*dFdvDict[var][iref]*parmDict[phfx+'Scale']*dervCor*ref[11]
    17891986                    if phfx+'Scale' in varylist:
    17901987                        dMdvh[varylist.index(phfx+'Scale')][iref] = w*ref[9]*dervCor
Note: See TracChangeset for help on using the changeset viewer.