Changeset 2944 for branch/2frame/GSASIImpsubs.py
 Timestamp:
 Jul 25, 2017 6:21:06 PM (4 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

branch/2frame/GSASIImpsubs.py
r2941 r2944 28 28 import numpy as np 29 29 import numpy.ma as ma 30 import numpy.linalg as nl 30 31 import GSASIIpath 31 32 GSASIIpath.SetVersionNumber("$Revision: 2895 $") … … 33 34 import GSASIIstrMath as G2stMth 34 35 36 sind = lambda x: np.sin(x*np.pi/180.) 35 37 cosd = lambda x: np.cos(x*np.pi/180.) 36 38 tand = lambda x: np.tan(x*np.pi/180.) 39 #asind = lambda x: 180.*np.arcsin(x)/np.pi 37 40 #acosd = lambda x: 180.*np.arccos(x)/np.pi 38 41 #atan2d = lambda y,x: 180.*np.arctan2(y,x)/np.pi … … 40 43 ncores = None 41 44 42 def InitMP( ):45 def InitMP(allowMP=True): 43 46 '''Called in routines to initialize use of Multiprocessing 44 47 ''' … … 46 49 if ncores is not None: return 47 50 useMP = False 51 if not allowMP: 52 print('Multiprocessing disabled') 53 ncores = 0 54 return 48 55 ncores = GSASIIpath.GetConfigValue('Multiprocessing_cores',1) 49 56 if ncores < 0: ncores = mp.cpu_count() … … 59 66 def InitDerivGlobals(im1,calcControls1,SGMT1,hfx1,phfx1,pfx1,G1,GB1,g1,SGData1, 60 67 parmDict1,wave1,shl1,x1,cw1,Ka21,A1,varylist1,dependentVars1, 61 dFdvDict1,lamRatio1,kRatio1,doPawley1 ):68 dFdvDict1,lamRatio1,kRatio1,doPawley1,pawleyLookup1): 62 69 '''Initialize for the computation of derivatives. Puts lots of junk into the global 63 70 namespace in this module, including the arrays for derivatives (when needed.) 64 71 ''' 65 72 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 67 74 im = im1 68 75 calcControls = calcControls1 … … 88 95 kRatio = kRatio1 89 96 doPawley = doPawley1 97 pawleyLookup = pawleyLookup1 90 98 # determine the parameters that will have derivatives computed only at end 91 99 global nonatomvarylist … … 164 172 dIdsh,dIdsp,dIdpola,dIdPO,dFdODF,dFdSA,dFdAb,dFdEx = G2stMth.GetIntensityDerv(refl,im,wave,Uniq,G,g,pfx,phfx,hfx,SGData,calcControls,parmDict) 165 173 pos = refl[5+im] 174 calcKa2 = False 166 175 if 'C' in calcControls[hfx+'histType']: 167 176 tanth = tand(pos/2.0) 168 177 costh = cosd(pos/2.0) 169 calcKa2 = False170 178 if Ka2: 171 179 pos2 = refl[5+im]+lamRatio*tanth # + 360/pi * Dlam/lam * tan(th) … … 222 230 names.update({hfx+'Absorption':[dFdAb,'int'],}) 223 231 else: #'T'OF 224 dpdA,dpdZ,dpdDC,dpdDA,dpdDB,dpdV = G etReflPosDerv(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) 225 233 names = {hfx+'Scale':[dIdsh,'int'],phfx+'Scale':[dIdsp,'int'], 226 234 hfx+'difC':[dpdDC,'pos'],hfx+'difA':[dpdDA,'pos'],hfx+'difB':[dpdDB,'pos'], … … 345 353 #if calcKa2: 346 354 # 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 360 def 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 376 def 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 390 def 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 404 def 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 iBegiFin: #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,pos2fmin)) 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 434 def 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 iBegiFin: #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 ################################################################################ 454 def 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 471 def 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 478 def 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
Note: See TracChangeset
for help on using the changeset viewer.