Changeset 2944


Ignore:
Timestamp:
Jul 25, 2017 6:21:06 PM (4 years ago)
Author:
toby
Message:

Finally: multiprocessing for refinement

Location:
branch/2frame
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • branch/2frame/GSASIIdataGUI.py

    r2942 r2944  
    47234723            G2frame.dataWindow.SetDataSize()
    47244724
    4725     which repaints the window. For routines [such as GSASII.UpdatePeakGrid()]
     4725    which repaints the window. For routines [such as GSASIIpwdGUI.UpdatePeakGrid()]
    47264726    that are called repeatedly to update the entire contents of dataWindow
    47274727    themselves, it is important to add calls to
  • branch/2frame/GSASIImpsubs.py

    r2941 r2944  
    2828import numpy as np
    2929import numpy.ma as ma
     30import numpy.linalg as nl
    3031import GSASIIpath
    3132GSASIIpath.SetVersionNumber("$Revision: 2895 $")
     
    3334import GSASIIstrMath as G2stMth
    3435
     36sind = lambda x: np.sin(x*np.pi/180.)
    3537cosd = lambda x: np.cos(x*np.pi/180.)
    3638tand = lambda x: np.tan(x*np.pi/180.)
     39#asind = lambda x: 180.*np.arcsin(x)/np.pi
    3740#acosd = lambda x: 180.*np.arccos(x)/np.pi
    3841#atan2d = lambda y,x: 180.*np.arctan2(y,x)/np.pi
     
    4043ncores = None
    4144
    42 def InitMP():
     45def InitMP(allowMP=True):
    4346    '''Called in routines to initialize use of Multiprocessing
    4447    '''
     
    4649    if ncores is not None: return
    4750    useMP = False
     51    if not allowMP:
     52        print('Multiprocessing disabled')
     53        ncores = 0
     54        return
    4855    ncores = GSASIIpath.GetConfigValue('Multiprocessing_cores',-1)
    4956    if ncores < 0: ncores = mp.cpu_count()
     
    5966def InitDerivGlobals(im1,calcControls1,SGMT1,hfx1,phfx1,pfx1,G1,GB1,g1,SGData1,
    6067                     parmDict1,wave1,shl1,x1,cw1,Ka21,A1,varylist1,dependentVars1,
    61                      dFdvDict1,lamRatio1,kRatio1,doPawley1):
     68                     dFdvDict1,lamRatio1,kRatio1,doPawley1,pawleyLookup1):
    6269    '''Initialize for the computation of derivatives. Puts lots of junk into the global
    6370    namespace in this module, including the arrays for derivatives (when needed.)
    6471    '''
    6572    global im,calcControls,SGMT,hfx,phfx,pfx,G,GB,g,SGData,parmDict,wave,shl,x,cw,Ka2,A
    66     global varylist,dependentVars,dFdvDict,lamRatio,kRatio,doPawley
     73    global varylist,dependentVars,dFdvDict,lamRatio,kRatio,doPawley,pawleyLookup
    6774    im = im1
    6875    calcControls = calcControls1
     
    8895    kRatio = kRatio1
    8996    doPawley = doPawley1
     97    pawleyLookup = pawleyLookup1
    9098    # determine the parameters that will have derivatives computed only at end
    9199    global nonatomvarylist
     
    164172    dIdsh,dIdsp,dIdpola,dIdPO,dFdODF,dFdSA,dFdAb,dFdEx = G2stMth.GetIntensityDerv(refl,im,wave,Uniq,G,g,pfx,phfx,hfx,SGData,calcControls,parmDict)
    165173    pos = refl[5+im]
     174    calcKa2 = False
    166175    if 'C' in calcControls[hfx+'histType']:
    167176        tanth = tand(pos/2.0)
    168177        costh = cosd(pos/2.0)
    169         calcKa2 = False
    170178        if Ka2:
    171179            pos2 = refl[5+im]+lamRatio*tanth       # + 360/pi * Dlam/lam * tan(th)
     
    222230            names.update({hfx+'Absorption':[dFdAb,'int'],})
    223231    else:   #'T'OF
    224         dpdA,dpdZ,dpdDC,dpdDA,dpdDB,dpdV = GetReflPosDerv(refl,im,0.0,A,pfx,hfx,calcControls,parmDict)
     232        dpdA,dpdZ,dpdDC,dpdDA,dpdDB,dpdV = G2stMth.GetReflPosDerv(refl,im,0.0,A,pfx,hfx,calcControls,parmDict)
    225233        names = {hfx+'Scale':[dIdsh,'int'],phfx+'Scale':[dIdsp,'int'],
    226234            hfx+'difC':[dpdDC,'pos'],hfx+'difA':[dpdDA,'pos'],hfx+'difB':[dpdDB,'pos'],
     
    345353           #if calcKa2:
    346354           #    depDerivDict[name][iBeg:iFin] += dFdvDict[name][iref]*corr2
     355
     356################################################################################
     357# Fobs Squared computation
     358################################################################################       
     359#x,ratio,shl,xB,xF,im,lamRatio,kRatio,xMask,Ka2
     360def InitFobsSqGlobals(x1,ratio1,shl1,xB1,xF1,im1,lamRatio1,kRatio1,xMask1,Ka21):
     361    '''Initialize for the computation of Fobs Squared for powder histograms.
     362    Puts lots of junk into the global namespace in this module.
     363    '''
     364    global x,ratio,shl,xB,xF,im,lamRatio,kRatio,xMask,Ka2
     365    x = ma.getdata(x1)
     366    ratio = ratio1
     367    shl = shl1
     368    xB = xB1
     369    xF = xF1
     370    im = im1
     371    lamRatio = lamRatio1
     372    kRatio = kRatio1
     373    xMask = xMask1
     374    Ka2 = Ka21
     375
     376def ComputeFobsSqCWbatch(profList):
     377    sInt = 0
     378    resList = []
     379    for refl,iref in profList:
     380        icod = ComputeFobsSqCW(refl,iref)
     381        if type(icod) is tuple:
     382            resList.append((icod[0],iref))
     383            sInt += icod[1]
     384        elif icod == -1:
     385            res.append((None,iref))
     386        elif icod == -2:
     387            break
     388    return sInt,resList
     389
     390def ComputeFobsSqTOFbatch(profList):
     391    sInt = 0
     392    resList = []
     393    for refl,iref in profList:
     394        icod = ComputeFobsSqTOF(refl,iref)
     395        if type(icod) is tuple:
     396            resList.append((icod[0],iref))
     397            sInt += icod[1]
     398        elif icod == -1:
     399            res.append((None,iref))
     400        elif icod == -2:
     401            break
     402    return sInt,resList
     403       
     404def ComputeFobsSqCW(refl,iref):
     405    yp = np.zeros(len(x)) # not masked
     406    sInt = 0
     407    refl8im = 0
     408    Wd,fmin,fmax = G2pwd.getWidthsCW(refl[5+im],refl[6+im],refl[7+im],shl)
     409    iBeg = max(xB,np.searchsorted(x,refl[5+im]-fmin))
     410    iFin = max(xB,min(np.searchsorted(x,refl[5+im]+fmax),xF))
     411    iFin2 = iFin
     412    if not iBeg+iFin:       #peak below low limit - skip peak
     413        return 0
     414    if ma.all(xMask[iBeg:iFin]):    #peak entirely masked - skip peak
     415        return -1
     416    elif not iBeg-iFin:     #peak above high limit - done
     417        return -2
     418    elif iBeg < iFin:
     419        yp[iBeg:iFin] = refl[11+im]*refl[9+im]*G2pwd.getFCJVoigt3(
     420            refl[5+im],refl[6+im],refl[7+im],shl,x[iBeg:iFin])
     421        sInt = refl[11+im]*refl[9+im]
     422        if Ka2:
     423            pos2 = refl[5+im]+lamRatio*tand(refl[5+im]/2.0)       # + 360/pi * Dlam/lam * tan(th)
     424            Wd,fmin,fmax = G2pwd.getWidthsCW(pos2,refl[6+im],refl[7+im],shl)
     425            iBeg2 = max(xB,np.searchsorted(x,pos2-fmin))
     426            iFin2 = min(np.searchsorted(x,pos2+fmax),xF)
     427            if iFin2 > iBeg2:
     428                yp[iBeg2:iFin2] += refl[11+im]*refl[9+im]*kRatio*G2pwd.getFCJVoigt3(
     429                    pos2,refl[6+im],refl[7+im],shl,x[iBeg2:iFin2])
     430                sInt *= 1.+kRatio
     431    refl8im = np.sum(np.where(ratio[iBeg:iFin2]>0.,yp[iBeg:iFin2]*ratio[iBeg:iFin2]/(refl[11+im]*(1.+kRatio)),0.0))
     432    return refl8im,sInt
     433
     434def ComputeFobsSqTOF(refl,iref):
     435    yp = np.zeros(len(x)) # not masked
     436    refl8im = 0
     437    Wd,fmin,fmax = G2pwd.getWidthsTOF(refl[5+im],refl[12+im],refl[13+im],refl[6+im],refl[7+im])
     438    iBeg = max(xB,np.searchsorted(x,refl[5+im]-fmin))
     439    iFin = max(xB,min(np.searchsorted(x,refl[5+im]+fmax),xF))
     440    if not iBeg+iFin:       #peak below low limit - skip peak
     441        return 0
     442    if ma.all(xMask[iBeg:iFin]):    #peak entirely masked - skip peak
     443        return -1
     444    elif not iBeg-iFin:     #peak above high limit - done
     445        return -2
     446    if iBeg < iFin:
     447        yp[iBeg:iFin] = refl[11+im]*refl[9+im]*G2pwd.getEpsVoigt(
     448            refl[5+im],refl[12+im],refl[13+im],refl[6+im],refl[7+im],x[iBeg:iFin])
     449    refl8im = np.sum(np.where(ratio[iBeg:iFin]>0.,yp[iBeg:iFin]*ratio[iBeg:iFin]/refl[11+im],0.0))
     450    return refl8im,refl[11+im]*refl[9+im]
     451################################################################################
     452# Powder Profile computation
     453################################################################################       
     454def InitPwdrProfGlobals(im1,shl1,x1):
     455    '''Initialize for the computation of Fobs Squared for powder histograms.
     456    Puts lots of junk into the global namespace in this module.
     457    '''
     458    global im,shl,x
     459    im = im1
     460    shl = shl1
     461    x = ma.getdata(x1)
     462    global cw
     463    cw = np.diff(x)
     464    cw = np.append(cw,cw[-1])
     465    # create local copies of ycalc array
     466    if useMP:
     467        global yc
     468        yc = np.zeros_like(x1)
     469
     470
     471def ComputePwdrProfCW(profList):
     472    'Compute the peaks profile for a set of CW peaks and add into the yc array'
     473    for pos,refl,iBeg,iFin,kRatio in profList:
     474        yc[iBeg:iFin] += refl[11+im]*refl[9+im]*kRatio*G2pwd.getFCJVoigt3(
     475            pos,refl[6+im],refl[7+im],shl,x[iBeg:iFin])
     476    return yc
     477
     478def ComputePwdrProfTOF(profList):
     479    'Compute the peaks profile for a set of TOF peaks and add into the yc array'
     480    for pos,refl,iBeg,iFin in profList:
     481        yc[iBeg:iFin] += refl[11+im]*refl[9+im]*G2pwd.getEpsVoigt(
     482            pos,refl[12+im],refl[13+im],refl[6+im],refl[7+im],x[iBeg:iFin])/cw[iBeg:iFin]
     483    return yc
     484   
  • branch/2frame/GSASIIstrMath.py

    r2935 r2944  
    1313import time
    1414import copy
     15import multiprocessing as mp
    1516import numpy as np
    1617import numpy.ma as ma
     
    2627import GSASIImath as G2mth
    2728import GSASIIobj as G2obj
     29import GSASIImpsubs as G2mp
    2830
    2931sind = lambda x: np.sin(x*np.pi/180.)
     
    29242926def GetFobsSq(Histograms,Phases,parmDict,calcControls):
    29252927    'Compute the observed structure factors for Powder histograms'
    2926     #starttime = time.time(); print 'start GetFobsSq'
     2928    #if GSASIIpath.GetConfigValue('debug'):
     2929    #    starttime = time.time(); print 'start GetFobsSq'
     2930    G2mp.InitMP()
    29272931    histoList = Histograms.keys()
    29282932    histoList.sort()
     2933    Ka2 = shl = lamRatio = kRatio = None
    29292934    for histogram in histoList:
    29302935        if 'PWDR' in histogram[:4]:
     
    29482953            ymb = np.where(ymb,ymb,1.0)
    29492954            ycmb = np.array(yc-yb)
    2950             ratio = 1./np.where(ycmb,ycmb/ymb,1.e10)         
     2955            ratio = 1./np.where(ycmb,ycmb/ymb,1.e10)
    29512956            refLists = Histogram['Reflection Lists']
    29522957            for phase in refLists:
     
    29602965                phfx = '%d:%d:'%(pId,hId)
    29612966                refDict = refLists[phase]
     2967                sumInt = 0.0
     2968                nExcl = 0
     2969                useMP = G2mp.useMP and len(refDict['RefList']) > 100
     2970                if useMP: # yes, create a set of initialized Python processes
     2971                    MPpool = mp.Pool(G2mp.ncores,G2mp.InitFobsSqGlobals,
     2972                                    [x,ratio,shl,xB,xF,im,lamRatio,kRatio,xMask,Ka2])
     2973                    profArgs = [[] for i in range(G2mp.ncores)]
     2974                else:
     2975                    G2mp.InitFobsSqGlobals(x,ratio,shl,xB,xF,im,lamRatio,kRatio,xMask,Ka2)
     2976                if 'C' in calcControls[hfx+'histType']:
     2977                    # are we multiprocessing?
     2978                    for iref,refl in enumerate(refDict['RefList']):
     2979                        if useMP:
     2980                            profArgs[iref%G2mp.ncores].append((refl,iref))
     2981                        else:
     2982                            icod= G2mp.ComputeFobsSqCW(refl,iref)
     2983                            if type(icod) is tuple:
     2984                                refl[8+im] = icod[0]
     2985                                sumInt += icod[1]
     2986                                if parmDict[phfx+'LeBail']: refl[9+im] = refl[8+im]
     2987                            elif icod == -1:
     2988                                refl[3+im] *= -1
     2989                                nExcl += 1
     2990                            elif icod == -2:
     2991                                break
     2992                    if useMP:
     2993                        for sInt,resList in MPpool.imap_unordered(G2mp.ComputeFobsSqCWbatch,profArgs):
     2994                            sumInt += sInt
     2995                            for refl8im,irefl in resList:
     2996                                if refl8im is None:
     2997                                    refDict['RefList'][irefl][3+im] *= -1
     2998                                    nExcl += 1
     2999                                else:
     3000                                    refDict['RefList'][irefl][8+im] = refl8im
     3001                                    if parmDict[phfx+'LeBail']:
     3002                                        refDict['RefList'][irefl][9+im] = refDict['RefList'][irefl][8+im]
     3003                        MPpool.terminate()
     3004                elif 'T' in calcControls[hfx+'histType']:
     3005                    for iref,refl in enumerate(refDict['RefList']):
     3006                        if useMP:
     3007                            profArgs[iref%G2mp.ncores].append((refl,iref))
     3008                        else:
     3009                            icod= G2mp.ComputeFobsSqTOF(refl,iref)
     3010                            if type(icod) is tuple:
     3011                                refl[8+im] = icod[0]
     3012                                sumInt += icod[1]
     3013                                if parmDict[phfx+'LeBail']: refl[9+im] = refl[8+im]
     3014                            elif icod == -1:
     3015                                refl[3+im] *= -1
     3016                                nExcl += 1
     3017                            elif icod == -2:
     3018                                break
     3019                    if useMP:
     3020                        for sInt,resList in MPpool.imap_unordered(G2mp.ComputeFobsSqTOFbatch,profArgs):
     3021                            sumInt += sInt
     3022                            for refl8im,irefl in resList:
     3023                                if refl8im is None:
     3024                                    refDict['RefList'][irefl][3+im] *= -1
     3025                                    nExcl += 1
     3026                                else:
     3027                                    refDict['RefList'][irefl][8+im] = refl8im
     3028                                    if parmDict[phfx+'LeBail']:
     3029                                        refDict['RefList'][irefl][9+im] = refDict['RefList'][irefl][8+im]
    29623030                sumFo = 0.0
    29633031                sumdF = 0.0
    29643032                sumFosq = 0.0
    29653033                sumdFsq = 0.0
    2966                 sumInt = 0.0
    2967                 nExcl = 0
    2968                 for refl in refDict['RefList']:
    2969                     if 'C' in calcControls[hfx+'histType']:
    2970                         yp = np.zeros_like(yb)
    2971                         Wd,fmin,fmax = G2pwd.getWidthsCW(refl[5+im],refl[6+im],refl[7+im],shl)
    2972                         iBeg = max(xB,np.searchsorted(x,refl[5+im]-fmin))
    2973                         iFin = max(xB,min(np.searchsorted(x,refl[5+im]+fmax),xF))
    2974                         iFin2 = iFin
    2975                         if not iBeg+iFin:       #peak below low limit - skip peak
    2976                             continue
    2977                         if ma.all(xMask[iBeg:iFin]):    #peak entirely masked - skip peak
    2978                             refl[3+im] *= -1
    2979                             nExcl += 1
    2980                             continue
    2981                         elif not iBeg-iFin:     #peak above high limit - done
    2982                             break
    2983                         elif iBeg < iFin:
    2984                             yp[iBeg:iFin] = refl[11+im]*refl[9+im]*G2pwd.getFCJVoigt3(refl[5+im],refl[6+im],refl[7+im],shl,ma.getdata(x[iBeg:iFin]))    #>90% of time spent here
    2985                             sumInt += refl[11+im]*refl[9+im]
    2986                             if Ka2:
    2987                                 pos2 = refl[5+im]+lamRatio*tand(refl[5+im]/2.0)       # + 360/pi * Dlam/lam * tan(th)
    2988                                 Wd,fmin,fmax = G2pwd.getWidthsCW(pos2,refl[6+im],refl[7+im],shl)
    2989                                 iBeg2 = max(xB,np.searchsorted(x,pos2-fmin))
    2990                                 iFin2 = min(np.searchsorted(x,pos2+fmax),xF)
    2991                                 if iFin2 > iBeg2:
    2992                                     yp[iBeg2:iFin2] += refl[11+im]*refl[9+im]*kRatio*G2pwd.getFCJVoigt3(pos2,refl[6+im],refl[7+im],shl,ma.getdata(x[iBeg2:iFin2]))        #and here
    2993                                     sumInt += refl[11+im]*refl[9+im]*kRatio
    2994                             refl[8+im] = np.sum(np.where(ratio[iBeg:iFin2]>0.,yp[iBeg:iFin2]*ratio[iBeg:iFin2]/(refl[11+im]*(1.+kRatio)),0.0))
    2995                             if parmDict[phfx+'LeBail']:
    2996                                 refl[9+im] = refl[8+im]
    2997                                
    2998                     elif 'T' in calcControls[hfx+'histType']:
    2999                         yp = np.zeros_like(yb)
    3000                         Wd,fmin,fmax = G2pwd.getWidthsTOF(refl[5+im],refl[12+im],refl[13+im],refl[6+im],refl[7+im])
    3001                         iBeg = max(xB,np.searchsorted(x,refl[5+im]-fmin))
    3002                         iFin = max(xB,min(np.searchsorted(x,refl[5+im]+fmax),xF))
    3003                         if not iBeg+iFin:       #peak below low limit - skip peak
    3004                             continue
    3005                         if ma.all(xMask[iBeg:iFin]):    #peak entirely masked - skip peak
    3006                             refl[3+im] *= -1
    3007                             nExcl += 1
    3008                             continue
    3009                         elif not iBeg-iFin:     #peak above high limit - done
    3010                             break
    3011                         if iBeg < iFin:
    3012                             yp[iBeg:iFin] = refl[11+im]*refl[9+im]*G2pwd.getEpsVoigt(refl[5+im],refl[12+im],refl[13+im],refl[6+im],refl[7+im],ma.getdata(x[iBeg:iFin]))  #>90% of time spent here
    3013                             refl[8+im] = np.sum(np.where(ratio[iBeg:iFin]>0.,yp[iBeg:iFin]*ratio[iBeg:iFin]/refl[11+im],0.0))
    3014                             if parmDict[phfx+'LeBail']:
    3015                                 refl[9+im] = refl[8+im]
    3016                             sumInt += refl[11+im]*refl[9+im]
     3034                for iref,refl in enumerate(refDict['RefList']):
    30173035                    Fo = np.sqrt(np.abs(refl[8+im]))
    30183036                    Fc = np.sqrt(np.abs(refl[9]+im))
     
    30333051            Histogram = Histograms[histogram]
    30343052            Histogram['Residuals']['hId'] = Histograms[histogram]['hId']
    3035     #print 'end GetFobsSq t=',time.time()-starttime
    3036                
     3053    #if GSASIIpath.GetConfigValue('debug'):
     3054    #    print 'end GetFobsSq t=',time.time()-starttime
     3055    #    print "Rf^2",[(i,Histogram['Residuals'][i]) for i in Histogram['Residuals'] if 'Rf^2' in i]
     3056
    30373057def getPowderProfile(parmDict,x,varylist,Histogram,Phases,calcControls,pawleyLookup):
    30383058    'Computes the powder pattern for a histogram based on contributions from all used phases'
    3039    
    3040     #starttime = time.time(); print 'start getPowderProfile'
    3041    
     3059    #if GSASIIpath.GetConfigValue('debug'):
     3060    #    starttime = time.time(); print 'start getPowderProfile'
     3061    G2mp.InitMP()   
     3062           
    30423063    def GetReflSigGamCW(refl,im,wave,G,GB,phfx,calcControls,parmDict):
    30433064        U = parmDict[hfx+'U']
     
    30673088        bet = parmDict[hfx+'beta-0']+parmDict[hfx+'beta-1']/refl[4+im]**4+parmDict[hfx+'beta-q']/refl[4+im]**2
    30683089        return alp,bet
    3069        
     3090   
    30703091    hId = Histogram['hId']
    30713092    hfx = ':%d:'%(hId)
     
    30753096    cw = np.diff(x)
    30763097    cw = np.append(cw,cw[-1])
    3077        
     3098   
     3099    shl = wave = kRatio = lamRatio = None # variables that are passed and must be initialized even when not needed
    30783100    if 'C' in calcControls[hfx+'histType']:   
    30793101        shl = max(parmDict[hfx+'SH/L'],0.002)
     
    31153137                StructureFactor2(refDict,G,hfx,pfx,SGData,calcControls,parmDict)
    31163138        badPeak = False
    3117         for iref,refl in enumerate(refDict['RefList']):
    3118             if 'C' in calcControls[hfx+'histType']:
     3139
     3140        # are we multiprocessing?
     3141        useMP = G2mp.useMP and len(refDict['RefList']) > 100
     3142        if useMP: # yes, create a set of initialized Python processes
     3143            MPpool = mp.Pool(G2mp.ncores,G2mp.InitPwdrProfGlobals,[im,shl,x])
     3144            profArgs = [[] for i in range(G2mp.ncores)]
     3145        if 'C' in calcControls[hfx+'histType']:
     3146            for iref,refl in enumerate(refDict['RefList']):
    31193147                if im:
    31203148                    h,k,l,m = refl[:4]
     
    31483176                    badPeak = True
    31493177                    continue
    3150                 yc[iBeg:iFin] += refl[11+im]*refl[9+im]*G2pwd.getFCJVoigt3(refl[5+im],refl[6+im],refl[7+im],shl,ma.getdata(x[iBeg:iFin]))    #>90% of time spent here
     3178                if useMP:
     3179                    profArgs[iref%G2mp.ncores].append((refl[5+im],refl,iBeg,iFin,1.))
     3180                else:
     3181                    yc[iBeg:iFin] += refl[11+im]*refl[9+im]*G2pwd.getFCJVoigt3(
     3182                        refl[5+im],refl[6+im],refl[7+im],shl,ma.getdata(x[iBeg:iFin]))
    31513183                if Ka2:
    31523184                    pos2 = refl[5+im]+lamRatio*tand(refl[5+im]/2.0)       # + 360/pi * Dlam/lam * tan(th)
     
    31603192                    elif iBeg > iFin:   #bad peak coeff - skip
    31613193                        continue
    3162                     yc[iBeg:iFin] += refl[11+im]*refl[9+im]*kRatio*G2pwd.getFCJVoigt3(pos2,refl[6+im],refl[7+im],shl,ma.getdata(x[iBeg:iFin]))        #and here
    3163             elif 'T' in calcControls[hfx+'histType']:
     3194                    if useMP:
     3195                        profArgs[iref%G2mp.ncores].append((pos2,refl,iBeg,iFin,kRatio))
     3196                    else:
     3197                        yc[iBeg:iFin] += refl[11+im]*refl[9+im]*kRatio*G2pwd.getFCJVoigt3(
     3198                            pos2,refl[6+im],refl[7+im],shl,ma.getdata(x[iBeg:iFin]))
     3199        elif 'T' in calcControls[hfx+'histType']:
     3200            for iref,refl in enumerate(refDict['RefList']):
    31643201                h,k,l = refl[:3]
    31653202                Uniq = np.inner(refl[:3],SGMT)
    3166                 refl[5+im] = GetReflPos(refl,im,0.0,A,pfx,hfx,calcControls,parmDict)         #corrected reflection position - #TODO - what about tabluated offset?
     3203                refl[5+im] = GetReflPos(refl,im,0.0,A,pfx,hfx,calcControls,parmDict)         #corrected reflection position
    31673204                Lorenz = sind(abs(parmDict[hfx+'2-theta'])/2)*refl[4+im]**4                                                #TOF Lorentz correction
    31683205#                refl[5+im] += GetHStrainShift(refl,im,SGData,phfx,hfx,calcControls,parmDict)               #apply hydrostatic strain shift
    31693206                refl[6+im:8+im] = GetReflSigGamTOF(refl,im,G,GB,phfx,calcControls,parmDict)    #peak sig & gam
    3170                 refl[12+im:14+im] = GetReflAlpBet(refl,im,hfx,parmDict)             #TODO - skip if alp, bet tabulated?
     3207                refl[12+im:14+im] = GetReflAlpBet(refl,im,hfx,parmDict)
    31713208                refl[11+im],refl[15+im],refl[16+im],refl[17+im] = GetIntensityCorr(refl,im,Uniq,G,g,pfx,phfx,hfx,SGData,calcControls,parmDict)
    31723209                refl[11+im] *= Vst*Lorenz
     
    31913228                    badPeak = True
    31923229                    continue
    3193                 yc[iBeg:iFin] += refl[11+im]*refl[9+im]*G2pwd.getEpsVoigt(refl[5+im],refl[12+im],refl[13+im],refl[6+im],refl[7+im],ma.getdata(x[iBeg:iFin]))/cw[iBeg:iFin]
    3194 #        print 'profile calc time: %.3fs'%(time.time()-time0)
     3230                if useMP:
     3231                    profArgs[iref%G2mp.ncores].append((refl[5+im],refl,iBeg,iFin))
     3232                else:
     3233                    yc[iBeg:iFin] += refl[11+im]*refl[9+im]*G2pwd.getEpsVoigt(refl[5+im],refl[12+im],refl[13+im],refl[6+im],refl[7+im],ma.getdata(x[iBeg:iFin]))/cw[iBeg:iFin]
     3234        if useMP and 'C' in calcControls[hfx+'histType']:
     3235            for y in MPpool.imap_unordered(G2mp.ComputePwdrProfCW,profArgs):
     3236                yc += y
     3237            MPpool.terminate()
     3238        elif useMP:
     3239            for y in MPpool.imap_unordered(G2mp.ComputePwdrProfTOF,profArgs):
     3240                yc += y
     3241            MPpool.terminate()
    31953242    if badPeak:
    31963243        print 'ouch #4 bad profile coefficients yield negative peak width; some reflections skipped'
    3197     #print 'end getPowderProfile t=',time.time()-starttime
     3244    #if GSASIIpath.GetConfigValue('debug'):
     3245    #    print 'end getPowderProfile t=',time.time()-starttime
    31983246    return yc,yb
    31993247   
     
    32023250    refined parameters
    32033251    '''
    3204     #starttime = time.time(); print 'start getPowderProfileDerv'
    3205    
    3206     def cellVaryDerv(pfx,SGData,dpdA):
    3207         if SGData['SGLaue'] in ['-1',]:
    3208             return [[pfx+'A0',dpdA[0]],[pfx+'A1',dpdA[1]],[pfx+'A2',dpdA[2]],
    3209                 [pfx+'A3',dpdA[3]],[pfx+'A4',dpdA[4]],[pfx+'A5',dpdA[5]]]
    3210         elif SGData['SGLaue'] in ['2/m',]:
    3211             if SGData['SGUniq'] == 'a':
    3212                 return [[pfx+'A0',dpdA[0]],[pfx+'A1',dpdA[1]],[pfx+'A2',dpdA[2]],[pfx+'A5',dpdA[5]]]
    3213             elif SGData['SGUniq'] == 'b':
    3214                 return [[pfx+'A0',dpdA[0]],[pfx+'A1',dpdA[1]],[pfx+'A2',dpdA[2]],[pfx+'A4',dpdA[4]]]
    3215             else:
    3216                 return [[pfx+'A0',dpdA[0]],[pfx+'A1',dpdA[1]],[pfx+'A2',dpdA[2]],[pfx+'A3',dpdA[3]]]
    3217         elif SGData['SGLaue'] in ['mmm',]:
    3218             return [[pfx+'A0',dpdA[0]],[pfx+'A1',dpdA[1]],[pfx+'A2',dpdA[2]]]
    3219         elif SGData['SGLaue'] in ['4/m','4/mmm']:
    3220             return [[pfx+'A0',dpdA[0]],[pfx+'A2',dpdA[2]]]
    3221         elif SGData['SGLaue'] in ['6/m','6/mmm','3m1', '31m', '3']:
    3222             return [[pfx+'A0',dpdA[0]],[pfx+'A2',dpdA[2]]]
    3223         elif SGData['SGLaue'] in ['3R', '3mR']:
    3224             return [[pfx+'A0',dpdA[0]+dpdA[1]+dpdA[2]],[pfx+'A3',dpdA[3]+dpdA[4]+dpdA[5]]]                       
    3225         elif SGData['SGLaue'] in ['m3m','m3']:
    3226             return [[pfx+'A0',dpdA[0]]]
    3227            
     3252    #if GSASIIpath.GetConfigValue('debug'):
     3253    #    starttime = time.time()
     3254    #    print 'starting getPowderProfileDerv'
     3255    G2mp.InitMP()   
    32283256    # create a list of dependent variables and set up a dictionary to hold their derivatives
    32293257    dependentVars = G2mv.GetDependentVars()
     
    32603288    cw = np.append(cw,cw[-1])
    32613289    Ka2 = False #also for TOF!
     3290    shl = wave = kRatio = lamRatio = None # variables that are passed and must be initialized even when not needed
    32623291    if 'C' in calcControls[hfx+'histType']:   
    32633292        shl = max(parmDict[hfx+'SH/L'],0.002)
     
    32693298        else:
    32703299            wave = parmDict[hfx+'Lam']
    3271     #print '#1 getPowderProfileDerv t=',time.time()-starttime
    32723300    for phase in Histogram['Reflection Lists']:
    32733301        refDict = Histogram['Reflection Lists'][phase]
     
    32993327#            print 'sf-derv time %.3fs'%(time.time()-time0)
    33003328            ApplyRBModelDervs(dFdvDict,parmDict,rigidbodyDict,Phase)
    3301         #print '#2 getPowderProfileDerv t=',time.time()-starttime
    3302         # determine the parameters that will have derivatives computed only at end
    3303         nonatomvarylist = []
    3304         for name in varylist:
    3305             if '::RBV;' not in name:
    3306                 try:
    3307                     aname = name.split(pfx)[1][:2]
    3308                     if aname not in ['Af','dA','AU','RB','AM','Xs','Xc','Ys','Yc','Zs','Zc',    \
    3309                         'Tm','Xm','Ym','Zm','U1','U2','U3']: continue # skip anything not an atom or rigid body param
    3310                 except IndexError:
    3311                     continue
    3312             nonatomvarylist.append(name)
    3313         nonatomdependentVars = []
    3314         for name in dependentVars:
    3315             if '::RBV;' not in name:
    3316                 try:
    3317                     aname = name.split(pfx)[1][:2]
    3318                     if aname not in ['Af','dA','AU','RB','AM','Xs','Xc','Ys','Yc','Zs','Zc',    \
    3319                         'Tm','Xm','Ym','Zm','U1','U2','U3']: continue # skip anything not an atom or rigid body param
    3320                 except IndexError:
    3321                     continue
    3322             nonatomdependentVars.append(name)
    3323         #timelist =  10*[0.0]
    3324         #timestart = 10*[0.0]
    3325         #==========================================================================================
    3326         #==========================================================================================
     3329
     3330        # are we multiprocessing?
     3331        useMP = G2mp.useMP and len(refDict['RefList']) > 100
     3332        if useMP: # yes, create a set of initialized Python processes
     3333            MPpool = mp.Pool(G2mp.ncores,G2mp.InitDerivGlobals,
     3334                            [im,calcControls,SGMT,hfx,phfx,pfx,G,GB,g,SGData,parmDict,
     3335                               wave,shl,x,cw,Ka2,A,varylist,dependentVars,dFdvDict,
     3336                               lamRatio,kRatio,Phase['General'].get('doPawley'),
     3337                               pawleyLookup
     3338                            ])
     3339            derivArgs = [[] for i in range(G2mp.ncores)]
     3340        else: # no, initialize the G2mp module
     3341            G2mp.InitDerivGlobals(im,calcControls,SGMT,hfx,phfx,pfx,G,GB,g,SGData,parmDict,
     3342                              wave,shl,x,cw,Ka2,A,varylist,dependentVars,dFdvDict,
     3343                              lamRatio,kRatio,Phase['General'].get('doPawley'),
     3344                              pawleyLookup)
     3345        # loop over all reflections
    33273346        for iref,refl in enumerate(refDict['RefList']):
    3328             #timestart[0] = time.time()
    3329             if im:
    3330                 h,k,l,m = refl[:4]
    3331             else:
    3332                 h,k,l = refl[:3]
    3333             Uniq = np.inner(refl[:3],SGMT)
    3334             if 'T' in calcControls[hfx+'histType']:
    3335                 wave = refl[14+im]
    3336             dIdsh,dIdsp,dIdpola,dIdPO,dFdODF,dFdSA,dFdAb,dFdEx = GetIntensityDerv(refl,im,wave,Uniq,G,g,pfx,phfx,hfx,SGData,calcControls,parmDict)
    33373347            if 'C' in calcControls[hfx+'histType']:        #CW powder
    33383348                Wd,fmin,fmax = G2pwd.getWidthsCW(refl[5+im],refl[6+im],refl[7+im],shl)
     
    33453355            elif not iBeg-iFin:     #peak above high limit - done
    33463356                break
    3347             pos = refl[5+im]
    3348             #itim=0;timelist[itim] += time.time()-timestart[itim]; timestart[itim+1] = time.time()
    3349             if 'C' in calcControls[hfx+'histType']:
    3350                 tanth = tand(pos/2.0)
    3351                 costh = cosd(pos/2.0)
    3352                 lenBF = iFin-iBeg
    3353                 dMdpk = np.zeros(shape=(6,lenBF))
    3354                 dMdipk = G2pwd.getdFCJVoigt3(refl[5+im],refl[6+im],refl[7+im],shl,ma.getdata(x[iBeg:iFin]))
    3355                 for i in range(5):
    3356                     dMdpk[i] += 100.*cw[iBeg:iFin]*refl[11+im]*refl[9+im]*dMdipk[i]
    3357                 dervDict = {'int':dMdpk[0],'pos':dMdpk[1],'sig':dMdpk[2],'gam':dMdpk[3],'shl':dMdpk[4],'L1/L2':np.zeros_like(dMdpk[0])}
    3358                 if Ka2:
    3359                     pos2 = refl[5+im]+lamRatio*tanth       # + 360/pi * Dlam/lam * tan(th)
    3360                     iBeg2 = np.searchsorted(x,pos2-fmin)
    3361                     iFin2 = np.searchsorted(x,pos2+fmax)
    3362                     if iBeg2-iFin2:
    3363                         lenBF2 = iFin2-iBeg2
    3364                         dMdpk2 = np.zeros(shape=(6,lenBF2))
    3365                         dMdipk2 = G2pwd.getdFCJVoigt3(pos2,refl[6+im],refl[7+im],shl,ma.getdata(x[iBeg2:iFin2]))
    3366                         for i in range(5):
    3367                             dMdpk2[i] = 100.*cw[iBeg2:iFin2]*refl[11+im]*refl[9+im]*kRatio*dMdipk2[i]
    3368                         dMdpk2[5] = 100.*cw[iBeg2:iFin2]*refl[11+im]*dMdipk2[0]
    3369                         dervDict2 = {'int':dMdpk2[0],'pos':dMdpk2[1],'sig':dMdpk2[2],'gam':dMdpk2[3],'shl':dMdpk2[4],'L1/L2':dMdpk2[5]*refl[9]}
    3370             else:   #'T'OF
    3371                 lenBF = iFin-iBeg
    3372                 if lenBF < 0:   #bad peak coeff
    3373                     break
    3374                 dMdpk = np.zeros(shape=(6,lenBF))
    3375                 dMdipk = G2pwd.getdEpsVoigt(refl[5+im],refl[12+im],refl[13+im],refl[6+im],refl[7+im],ma.getdata(x[iBeg:iFin]))
    3376                 for i in range(6):
    3377                     dMdpk[i] += refl[11+im]*refl[9+im]*dMdipk[i]      #cw[iBeg:iFin]*
    3378                 dervDict = {'int':dMdpk[0],'pos':dMdpk[1],'alp':dMdpk[2],'bet':dMdpk[3],'sig':dMdpk[4],'gam':dMdpk[5]}           
    3379             #itim=1;timelist[itim] += time.time()-timestart[itim]; timestart[itim+1] = time.time()
    3380             if Phase['General'].get('doPawley'):
    3381                 dMdpw = np.zeros(len(x))
    3382                 try:
    3383                     if im:
    3384                         pIdx = pfx+'PWLref:'+str(pawleyLookup[pfx+'%d,%d,%d,%d'%(h,k,l,m)])
    3385                     else:
    3386                         pIdx = pfx+'PWLref:'+str(pawleyLookup[pfx+'%d,%d,%d'%(h,k,l)])
    3387                     idx = varylist.index(pIdx)
    3388                     dMdpw[iBeg:iFin] = dervDict['int']/refl[9+im]
    3389                     if Ka2: #not for TOF either
    3390                         dMdpw[iBeg2:iFin2] += dervDict2['int']/refl[9+im]
    3391                     dMdv[idx] = dMdpw
    3392                 except: # ValueError:
    3393                     pass
    3394             if 'C' in calcControls[hfx+'histType']:
    3395                 dpdA,dpdw,dpdZ,dpdSh,dpdTr,dpdX,dpdY,dpdV = GetReflPosDerv(refl,im,wave,A,pfx,hfx,calcControls,parmDict)
    3396                 names = {hfx+'Scale':[dIdsh,'int'],hfx+'Polariz.':[dIdpola,'int'],phfx+'Scale':[dIdsp,'int'],
    3397                     hfx+'U':[tanth**2,'sig'],hfx+'V':[tanth,'sig'],hfx+'W':[1.0,'sig'],
    3398                     hfx+'X':[1.0/costh,'gam'],hfx+'Y':[tanth,'gam'],hfx+'SH/L':[1.0,'shl'],
    3399                     hfx+'I(L2)/I(L1)':[1.0,'L1/L2'],hfx+'Zero':[dpdZ,'pos'],hfx+'Lam':[dpdw,'pos'],
    3400                     hfx+'Shift':[dpdSh,'pos'],hfx+'Transparency':[dpdTr,'pos'],hfx+'DisplaceX':[dpdX,'pos'],
    3401                     hfx+'DisplaceY':[dpdY,'pos'],}
    3402                 if 'Bragg' in calcControls[hfx+'instType']:
    3403                     names.update({hfx+'SurfRoughA':[dFdAb[0],'int'],
    3404                         hfx+'SurfRoughB':[dFdAb[1],'int'],})
    3405                 else:
    3406                     names.update({hfx+'Absorption':[dFdAb,'int'],})
    3407             else:   #'T'OF
    3408                 dpdA,dpdZ,dpdDC,dpdDA,dpdDB,dpdV = GetReflPosDerv(refl,im,0.0,A,pfx,hfx,calcControls,parmDict)
    3409                 names = {hfx+'Scale':[dIdsh,'int'],phfx+'Scale':[dIdsp,'int'],
    3410                     hfx+'difC':[dpdDC,'pos'],hfx+'difA':[dpdDA,'pos'],hfx+'difB':[dpdDB,'pos'],
    3411                     hfx+'Zero':[dpdZ,'pos'],hfx+'X':[refl[4+im],'gam'],hfx+'Y':[refl[4+im]**2,'gam'],
    3412                     hfx+'alpha':[1./refl[4+im],'alp'],hfx+'beta-0':[1.0,'bet'],hfx+'beta-1':[1./refl[4+im]**4,'bet'],
    3413                     hfx+'beta-q':[1./refl[4+im]**2,'bet'],hfx+'sig-0':[1.0,'sig'],hfx+'sig-1':[refl[4+im]**2,'sig'],
    3414                     hfx+'sig-2':[refl[4+im]**4,'sig'],hfx+'sig-q':[1./refl[4+im]**2,'sig'],
    3415                     hfx+'Absorption':[dFdAb,'int'],phfx+'Extinction':[dFdEx,'int'],}
    3416             #itim=2;timelist[itim] += time.time()-timestart[itim]; timestart[itim+1] = time.time()
    3417             for name in names:
    3418                 item = names[name]
    3419                 if name in varylist:
    3420                     dMdv[varylist.index(name)][iBeg:iFin] += item[0]*dervDict[item[1]]
    3421                     if Ka2 and iFin2-iBeg2:
    3422                         dMdv[varylist.index(name)][iBeg2:iFin2] += item[0]*dervDict2[item[1]]
    3423                 elif name in dependentVars:
    3424                     depDerivDict[name][iBeg:iFin] += item[0]*dervDict[item[1]]
    3425                     if Ka2 and iFin2-iBeg2:
    3426                         depDerivDict[name][iBeg2:iFin2] += item[0]*dervDict2[item[1]]
    3427             for iPO in dIdPO:
    3428                 if iPO in varylist:
    3429                     dMdv[varylist.index(iPO)][iBeg:iFin] += dIdPO[iPO]*dervDict['int']
    3430                     if Ka2 and iFin2-iBeg2:
    3431                         dMdv[varylist.index(iPO)][iBeg2:iFin2] += dIdPO[iPO]*dervDict2['int']
    3432                 elif iPO in dependentVars:
    3433                     depDerivDict[iPO][iBeg:iFin] += dIdPO[iPO]*dervDict['int']
    3434                     if Ka2 and iFin2-iBeg2:
    3435                         depDerivDict[iPO][iBeg2:iFin2] += dIdPO[iPO]*dervDict2['int']
    3436             #itim=3;timelist[itim] += time.time()-timestart[itim]; timestart[itim+1] = time.time()
    3437             for i,name in enumerate(['omega','chi','phi']):
    3438                 aname = pfx+'SH '+name
    3439                 if aname in varylist:
    3440                     dMdv[varylist.index(aname)][iBeg:iFin] += dFdSA[i]*dervDict['int']
    3441                     if Ka2 and iFin2-iBeg2:
    3442                         dMdv[varylist.index(aname)][iBeg2:iFin2] += dFdSA[i]*dervDict2['int']
    3443                 elif aname in dependentVars:
    3444                     depDerivDict[aname][iBeg:iFin] += dFdSA[i]*dervDict['int']
    3445                     if Ka2 and iFin2-iBeg2:
    3446                         depDerivDict[aname][iBeg2:iFin2] += dFdSA[i]*dervDict2['int']
    3447             for iSH in dFdODF:
    3448                 if iSH in varylist:
    3449                     dMdv[varylist.index(iSH)][iBeg:iFin] += dFdODF[iSH]*dervDict['int']
    3450                     if Ka2 and iFin2-iBeg2:
    3451                         dMdv[varylist.index(iSH)][iBeg2:iFin2] += dFdODF[iSH]*dervDict2['int']
    3452                 elif iSH in dependentVars:
    3453                     depDerivDict[iSH][iBeg:iFin] += dFdODF[iSH]*dervDict['int']
    3454                     if Ka2 and iFin2-iBeg2:
    3455                         depDerivDict[iSH][iBeg2:iFin2] += dFdODF[iSH]*dervDict2['int']
    3456             cellDervNames = cellVaryDerv(pfx,SGData,dpdA)
    3457             for name,dpdA in cellDervNames:
    3458                 if name in varylist:
    3459                     dMdv[varylist.index(name)][iBeg:iFin] += dpdA*dervDict['pos']
    3460                     if Ka2 and iFin2-iBeg2:
    3461                         dMdv[varylist.index(name)][iBeg2:iFin2] += dpdA*dervDict2['pos']
    3462                 elif name in dependentVars: #need to scale for mixed phase constraints?
    3463                     depDerivDict[name][iBeg:iFin] += dpdA*dervDict['pos']
    3464                     if Ka2 and iFin2-iBeg2:
    3465                         depDerivDict[name][iBeg2:iFin2] += dpdA*dervDict2['pos']
    3466             dDijDict = GetHStrainShiftDerv(refl,im,SGData,phfx,hfx,calcControls,parmDict)
    3467             for name in dDijDict:
    3468                 if name in varylist:
    3469                     dMdv[varylist.index(name)][iBeg:iFin] += dDijDict[name]*dervDict['pos']
    3470                     if Ka2 and iFin2-iBeg2:
    3471                         dMdv[varylist.index(name)][iBeg2:iFin2] += dDijDict[name]*dervDict2['pos']
    3472                 elif name in dependentVars:
    3473                     depDerivDict[name][iBeg:iFin] += dDijDict[name]*dervDict['pos']
    3474                     if Ka2 and iFin2-iBeg2:
    3475                         depDerivDict[name][iBeg2:iFin2] += dDijDict[name]*dervDict2['pos']
    3476             #itim=4;timelist[itim] += time.time()-timestart[itim]; timestart[itim+1] = time.time()
    3477             for i,name in enumerate([pfx+'mV0',pfx+'mV1',pfx+'mV2']):
    3478                 if name in varylist:
    3479                     dMdv[varylist.index(name)][iBeg:iFin] += dpdV[i]*dervDict['pos']
    3480                     if Ka2 and iFin2-iBeg2:
    3481                         dMdv[varylist.index(name)][iBeg2:iFin2] += dpdV[i]*dervDict2['pos']
    3482                 elif name in dependentVars:
    3483                     depDerivDict[name][iBeg:iFin] += dpdV[i]*dervDict['pos']
    3484                     if Ka2 and iFin2-iBeg2:
    3485                         depDerivDict[name][iBeg2:iFin2] += dpdV[i]*dervDict2['pos']
    3486             if 'C' in calcControls[hfx+'histType']:
    3487                 sigDict,gamDict = GetSampleSigGamDerv(refl,im,wave,G,GB,SGData,hfx,phfx,calcControls,parmDict)
    3488             else:   #'T'OF
    3489                 sigDict,gamDict = GetSampleSigGamDerv(refl,im,0.0,G,GB,SGData,hfx,phfx,calcControls,parmDict)
    3490             for name in gamDict:
    3491                 if name in varylist:
    3492                     dMdv[varylist.index(name)][iBeg:iFin] += gamDict[name]*dervDict['gam']
    3493                     if Ka2 and iFin2-iBeg2:
    3494                         dMdv[varylist.index(name)][iBeg2:iFin2] += gamDict[name]*dervDict2['gam']
    3495                 elif name in dependentVars:
    3496                     depDerivDict[name][iBeg:iFin] += gamDict[name]*dervDict['gam']
    3497                     if Ka2 and iFin2-iBeg2:
    3498                         depDerivDict[name][iBeg2:iFin2] += gamDict[name]*dervDict2['gam']
    3499             for name in sigDict:
    3500                 if name in varylist:
    3501                     dMdv[varylist.index(name)][iBeg:iFin] += sigDict[name]*dervDict['sig']
    3502                     if Ka2 and iFin2-iBeg2:
    3503                         dMdv[varylist.index(name)][iBeg2:iFin2] += sigDict[name]*dervDict2['sig']
    3504                 elif name in dependentVars:
    3505                     depDerivDict[name][iBeg:iFin] += sigDict[name]*dervDict['sig']
    3506                     if Ka2 and iFin2-iBeg2:
    3507                         depDerivDict[name][iBeg2:iFin2] += sigDict[name]*dervDict2['sig']
    3508             for name in ['BabA','BabU']:
    3509                 if refl[9+im]:
    3510                     if phfx+name in varylist:
    3511                         dMdv[varylist.index(phfx+name)][iBeg:iFin] += parmDict[phfx+'Scale']*dFdvDict[phfx+name][iref]*dervDict['int']/refl[9+im]
    3512                         if Ka2 and iFin2-iBeg2:
    3513                             dMdv[varylist.index(phfx+name)][iBeg2:iFin2] += parmDict[phfx+'Scale']*dFdvDict[phfx+name][iref]*dervDict2['int']/refl[9+im]
    3514                     elif phfx+name in dependentVars:                   
    3515                         depDerivDict[phfx+name][iBeg:iFin] += parmDict[phfx+'Scale']*dFdvDict[phfx+name][iref]*dervDict['int']/refl[9+im]
    3516                         if Ka2 and iFin2-iBeg2:
    3517                             depDerivDict[phfx+name][iBeg2:iFin2] += parmDict[phfx+'Scale']*dFdvDict[phfx+name][iref]*dervDict2['int']/refl[9+im]                 
    3518             #itim=5;timelist[itim] += time.time()-timestart[itim]; timestart[itim+1] = time.time()
    3519             if not Phase['General'].get('doPawley') and not parmDict[phfx+'LeBail']:
    3520                 #do atom derivatives -  for RB,F,X & U so far - how do I scale mixed phase constraints?
    3521                 corr = 0.
    3522                 corr2 = 0.
    3523                 if refl[9+im]:             
    3524                     corr = dervDict['int']/refl[9+im]
    3525                     #if Ka2 and iFin2-iBeg2:
    3526                     #    corr2 = dervDict2['int']/refl[9+im]
    3527                 #itim=6;timelist[itim] += time.time()-timestart[itim]; timestart[itim+1] = time.time()
    3528                 for name in nonatomvarylist:
    3529                     dMdv[varylist.index(name)][iBeg:iFin] += dFdvDict[name][iref]*corr
    3530                     if Ka2 and iFin2-iBeg2:
    3531                        dMdv[varylist.index(name)][iBeg2:iFin2] += dFdvDict[name][iref]*corr2
    3532                 #itim=7;timelist[itim] += time.time()-timestart[itim]; timestart[itim+1] = time.time()
    3533                 for name in nonatomdependentVars:
    3534                    depDerivDict[name][iBeg:iFin] += dFdvDict[name][iref]*corr
    3535                    if Ka2 and iFin2-iBeg2:
    3536                        depDerivDict[name][iBeg2:iFin2] += dFdvDict[name][iref]*corr2
    3537                 #itim=8;timelist[itim] += time.time()-timestart[itim]   
    3538     #        print 'profile derv time: %.3fs'%(time.time()-time0)
     3357            if useMP:
     3358                derivArgs[iref%G2mp.ncores].append((refl,iref,fmin,fmax,iBeg,iFin))
     3359            else:
     3360                result = G2mp.ComputeDeriv(refl,iref,fmin,fmax,iBeg,iFin,dMdv,depDerivDict)
     3361                if result: break
     3362        if useMP:
     3363            for dm,depd in MPpool.imap_unordered(G2mp.ComputeDerivMPbatch,derivArgs):
     3364                dMdv += dm
     3365                for i in depDerivDict: depDerivDict[i] += depd[i]
     3366            MPpool.terminate()
    35393367    # now process derivatives in constraints
    3540     #print '#3 getPowderProfileDerv t=',time.time()-starttime
    3541     #print timelist,sum(timelist)
    35423368    dMdv = ma.array(dMdv,mask=np.outer(np.ones(len(varylist)),ma.getmaskarray(x)))      #x is a MaskedArray!
    35433369    G2mv.Dict2Deriv(varylist,depDerivDict,dMdv)
    3544     #print 'end getPowderProfileDerv t=',time.time()-starttime
     3370    #if GSASIIpath.GetConfigValue('debug'):
     3371    #    print 'end getPowderProfileDerv t=',time.time()-starttime
    35453372    return dMdv
    35463373   
     
    40183845        M = np.concatenate((M,np.sqrt(pWt)*pVals))
    40193846    return M
    4020                        
     3847
Note: See TracChangeset for help on using the changeset viewer.