Changeset 3104 for trunk/GSASIIimage.py
 Timestamp:
 Sep 29, 2017 12:25:46 PM (5 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

trunk/GSASIIimage.py
r2836 r3104 50 50 npatand = lambda x: 180.*np.arctan(x)/np.pi 51 51 npatan2d = lambda y,x: 180.*np.arctan2(y,x)/np.pi 52 nxs = np.newaxis 52 53 debug = False 53 54 … … 813 814 return True 814 815 815 def Make2ThetaAzimuthMap(data,masks,iLim,jLim,t imes): #most expensive part of integration!816 def Make2ThetaAzimuthMap(data,masks,iLim,jLim,tamp,times): #most expensive part of integration! 816 817 'Needs a doc string' 817 818 #transforms 2D image from x,y space to 2theta,azimuth space based on detector orientation … … 821 822 822 823 tay,tax = np.mgrid[iLim[0]+0.5:iLim[1]+.5,jLim[0]+.5:jLim[1]+.5] #bin centers not corners 823 tax = np.asfarray(tax*scalex,dtype=np.float32) 824 tay = np.asfarray(tay*scaley,dtype=np.float32) 824 tax = np.asfarray(tax*scalex,dtype=np.float32).flatten() 825 tay = np.asfarray(tay*scaley,dtype=np.float32).flatten() 825 826 nI = iLim[1]iLim[0] 826 827 nJ = jLim[1]jLim[0] … … 828 829 #make position masks here 829 830 frame = masks['Frames'] 830 tam = ma.make_mask_none((nI ,nJ))831 tam = ma.make_mask_none((nI*nJ)) 831 832 if frame: 832 tamp = ma.make_mask_none((1024*1024)) 833 tamp = ma.make_mask(pm.polymask(nI*nJ,tax.flatten(), 834 tay.flatten(),len(frame),frame,tamp)[:nI*nJ])True #switch to exclude around frame 835 tam = ma.mask_or(tam.flatten(),tamp) 833 tam = ma.mask_or(tam,ma.make_mask(pm.polymask(nI*nJ,tax, 834 tay,len(frame),frame,tamp)[:nI*nJ])True) 836 835 polygons = masks['Polygons'] 837 836 for polygon in polygons: 838 837 if polygon: 839 tam p = ma.make_mask_none((1024*1024))840 tamp = ma.make_mask(pm.polymask(nI*nJ,tax.flatten(),841 tay.flatten(),len(polygon),polygon,tamp)[:nI*nJ])842 tam = ma.mask_or(tam.flatten(),tamp)838 tam = ma.mask_or(tam,ma.make_mask(pm.polymask(nI*nJ,tax, 839 tay,len(polygon),polygon,tamp)[:nI*nJ])) 840 for X,Y,rsq in masks['Points'].T: 841 tam = ma.mask_or(tam,ma.getmask(ma.masked_less((taxX)**2+(tayY)**2,rsq))) 843 842 if tam.shape: tam = np.reshape(tam,(nI,nJ)) 844 spots = masks['Points']845 for X,Y,diam in spots:846 tamp = ma.getmask(ma.masked_less((taxX)**2+(tayY)**2,(diam/2.)**2))847 tam = ma.mask_or(tam,tamp)848 843 times[0] += time.time()t0 849 844 t0 = time.time() 850 TA = np.array(GetTthAzmG( tax,tay,data)) #includes geom. corr. as dist**2/d0**2  most expensive step845 TA = np.array(GetTthAzmG(np.reshape(tax,(nI,nJ)),np.reshape(tay,(nI,nJ)),data)) #includes geom. corr. as dist**2/d0**2  most expensive step 851 846 times[1] += time.time()t0 852 847 TA[1] = np.where(TA[1]<0,TA[1]+360,TA[1]) … … 875 870 tabs = ma.compressed(ma.array(tabs.flatten(),mask=tam)) #ones  later used for absorption corr. 876 871 return tax,tay,taz,tad,tabs 877 872 878 873 def ImageIntegrate(image,data,masks,blkSize=128,returnN=False): 879 874 'Integrate an image; called from OnIntegrateAll and OnIntegrate in G2imgGUI' #for q, log(q) bins need data['binType'] … … 898 893 if 'SASD' in data['type']: 899 894 muT = np.log(muT)/2. #Transmission to 1/2 thickness muT 895 Masks = copy.deepcopy(masks) 896 Masks['Points'] = np.array(Masks['Points']).T #get spots as X,Y,R arrays 897 if np.any(masks['Points']): 898 Masks['Points'][2] = np.square(Masks['Points'][2]/2.) 900 899 NST = np.zeros(shape=(numAzms,numChans),order='F',dtype=np.float32) 901 900 H0 = np.zeros(shape=(numAzms,numChans),order='F',dtype=np.float32) … … 905 904 tbeg = time.time() 906 905 times = [0,0,0,0,0] 906 tamp = ma.make_mask_none((1024*1024)) #NB: this array size used in the fortran histogram2d 907 907 for iBlk in range(nYBlks): 908 908 iBeg = iBlk*blkSize … … 912 912 jFin = min(jBeg+blkSize,Nx) 913 913 # next is most expensive step! 914 TA,tam = Make2ThetaAzimuthMap(data, masks,(iBeg,iFin),(jBeg,jFin),times) #2theta & azimuth arrays & create position mask914 TA,tam = Make2ThetaAzimuthMap(data,Masks,(iBeg,iFin),(jBeg,jFin),tamp,times) #2theta & azimuth arrays & create position mask 915 915 Block = image[iBeg:iFin,jBeg:jFin] 916 916 t0 = time.time() 917 tax,tay,taz,tad,tabs = Fill2ThetaAzimuthMap( masks,TA,tam,Block) #and apply masks917 tax,tay,taz,tad,tabs = Fill2ThetaAzimuthMap(Masks,TA,tam,Block) #and apply masks 918 918 del TA; del tam 919 919 times[2] += time.time()t0
Note: See TracChangeset
for help on using the changeset viewer.