From af9ecfc181c12a6f1c198d6a469add0f909d8c64 Mon Sep 17 00:00:00 2001 From: Marco Cammarata Date: Fri, 13 Jan 2017 14:49:48 +0100 Subject: [PATCH] new submodules (peaks,cell) and new functions (like backgorundSubtraction (in azav.py) --- xray/__init__.py | 1 + xray/azav.py | 38 ++++++++++++++++------- xray/cell.py | 68 +++++++++++++++++++++++++++++++++++++++++ xray/dataReduction.py | 10 ++++--- xray/dspace.py | 70 +++++++++++++++++++++++++++++++++++++++++++ xray/id9.py | 27 ++++++++++------- xray/peaks.py | 20 +++++++++++++ xray/storage.py | 54 ++++++++++++++++++++++++++------- xray/utils.py | 27 ++++++++++++++++- 9 files changed, 278 insertions(+), 37 deletions(-) create mode 100644 xray/cell.py create mode 100644 xray/dspace.py create mode 100644 xray/peaks.py diff --git a/xray/__init__.py b/xray/__init__.py index 56e5b3c..3fff056 100644 --- a/xray/__init__.py +++ b/xray/__init__.py @@ -1,6 +1,7 @@ from . import storage from . import utils from . import mask +from . import peaks try: from . import azav except ImportError: diff --git a/xray/azav.py b/xray/azav.py index a28039a..2d96448 100644 --- a/xray/azav.py +++ b/xray/azav.py @@ -52,7 +52,14 @@ def ai_as_dict(ai): def ai_as_str(ai): """ ai is a pyFAI azimuthal intagrator""" - return "#" + str(ai).replace("\n","\n#") + s=[ "# Detector : %s" % ai.detector.name, + "# Pixel [um] : %.2fx%.2f" % (ai.pixel1*1e6,ai.pixel2*1e6), + "# Distance [mm] : %.3f" % (ai.dist*1e3), + "# Center [mm] : %.3f,%.3f" % (ai.poni1*1e3,ai.poni2*1e3), + "# Center [px] : %.3f,%.3f" % (ai.poni1/ai.pixel1,ai.poni2/ai.pixel2), + "# Wavelength [A] : %.5f" % (ai.wavelength*1e10), + "# rot[1,2,3] [rad]: %.3f,%.3f,%.3f" % (ai.rot1,ai.rot2,ai.rot3) ] + return "\n".join(s) def do1d(ai, imgs, mask = None, npt_radial = 600, method = 'csr',safe=True,dark=10., polCorr = 1): """ ai is a pyFAI azimuthal intagrator @@ -218,7 +225,9 @@ def leastsq_circle(x,y): residu = np.sum((Ri - R)**2) return xc, yc, R -def find_center(img,psize=100e-6,dist=0.1,wavelength=0.8e-10,**kwargs): +def find_center(img,psize=100e-6,dist=0.1,wavelength=0.8e-10,center=None,reference=None,**kwargs): + """ center is the initial centr (can be None) + reference is a reference position to be plot in 2D plots """ plt.ion() kw = dict( pixel1 = psize, pixel2 = psize, dist = dist,wavelength=wavelength ) kw.update(kwargs) @@ -232,23 +241,30 @@ def find_center(img,psize=100e-6,dist=0.1,wavelength=0.8e-10,**kwargs): ans = "" print("Enter 'end' when done") while ans != "end": - if ans == "": + if center is None: print("Click on beam center:") plt.sca(ax_img); # set figure to use for mouse interaction - xc,yc = plt.ginput()[0] - else: - xc,yc = map(float,ans.split(",")) - print("Selected center:",xc,yc) - ai.set_poni1(xc*psize) - ai.set_poni2(yc*psize) + center = plt.ginput()[0] + print("Selected center:",center) + ai.set_poni1(center[0]*psize) + ai.set_poni2(center[1]*psize) q,az,i = do2d(ai,img) mesh = ax_pyfai.pcolormesh(q,az,i) mesh.set_clim( *np.percentile(i,(2,95) ) ) - ax_pyfai.set_title(str( (xc,yc) )) + ax_pyfai.set_title(str(center)) + if reference is not None: ax_pyfai.axvline(reference) plt.pause(0.01) plt.draw() ans=input("Enter to continue with clinking or enter xc,yc values ") - print("Final values: (in pixels) %.3f %.3f"%(xc,yc)) + if ans == '': + center = None + else: + try: + center = list(map(float,ans.split(","))) + except Exception as e: + center = None + if center == []: center = None + print("Final values: (in pixels) %.3f %.3f"%(center[0],center[1])) return ai def average(fileOrFolder,delays=slice(None),scale=1,norm=None,returnAll=False,plot=False, diff --git a/xray/cell.py b/xray/cell.py new file mode 100644 index 0000000..cd28428 --- /dev/null +++ b/xray/cell.py @@ -0,0 +1,68 @@ +from __future__ import division,print_function +import numpy as np +from numpy import sin,cos + +class Triclinic(object): + def __init__(self,a=1,b=1,c=1,alpha=90,beta=90,gamma=90): + self.a = a + self.b = b + self.c = c + alpha = alpha*np.pi/180 + beta = beta*np.pi/180 + gamma = gamma*np.pi/180 + self.alpha = alpha + self.beta = beta + self.gamma = gamma + + self._s11 = b**2 * c**2 * sin(alpha)**2 + self._s22 = a**2 * c**2 * sin(beta)**2 + self._s33 = a**2 * b**2 * sin(gamma)**2 + self._s12 = a*b*c**2*(cos(alpha) * cos(beta) - cos(gamma)) + self._s23 = a**2*b*c*(cos(beta) * cos(gamma) - cos(alpha)) + self._s13 = a*b**2*c*(cos(gamma) * cos(alpha) - cos(beta)) + self.V = (a*b*c)*np.sqrt(1-cos(alpha)**2 - cos(beta)**2 - cos(gamma)**2 + 2*cos(alpha)*cos(beta)*cos(gamma)) + + def __call__(self,h,k,l): return self.q(h,k,l) + + def d(self,h,k,l): + temp = self._s11*h**2 + \ + self._s22*k**2 + \ + self._s33*l**2 + \ + 2*self._s12*h*k+ \ + 2*self._s23*k*l+ \ + 2*self._s13*h*l + d = self.V/np.sqrt(temp) + return d + + def q(self,h,k,l): + return 2*np.pi/self.d(h,k,l) + +class Orthorombic(Triclinic): + def __init__(self,a=1,b=1,c=1): + Triclinic.__init__(self,a=a,b=b,c=c,alpha=90,beta=90,gamma=90) + +class Monoclinic(object): + def __init__(self,a=1,b=1,c=1,beta=90.): + self.a = a + self.b = b + self.c = c + beta = beta/np.pi*180 + self.beta = beta + + self.V = (a*b*c) + + def __call__(self,h,k,l): return self.Q(h,k,l) + + def Q(self,h,k,l): + temp = h**2/self.a**2 + (k*sin(self.beta))**2/self.b**2+l**2/self.c**2+2*h*l*cos(self.beta)/self.a/self.c + d = 1/np.sqrt(temp) + print(d) + return 2*np.pi/d + + + +ti3o5_lambda = Triclinic(a = 9.83776, b = 3.78674, c = 9.97069, beta = 91.2567) +ti3o5_beta = Triclinic(a = 9.7382 , b = 3.8005 , c = 9.4333 , beta = 91.496) +#ti3o5_beta = Monoclinic(a = 9.7382 , b = 3.8005 , c = 9.4333 , beta = 91.496) +ti3o5_alpha = Triclinic(a = 9.8372, b = 3.7921, c = 9.9717) +#ti3o5_alpha1 = Orthorombic(a = 9.8372, b = 3.7921, c = 9.9717) diff --git a/xray/dataReduction.py b/xray/dataReduction.py index c960195..005605d 100644 --- a/xray/dataReduction.py +++ b/xray/dataReduction.py @@ -73,8 +73,10 @@ def averageScanPoints(scan,data,isRef=None,lpower=None,useRatio=False,\ if isRef.sum()>0: # create a copy (subtractReferences works in place) diff = subtractReferences(data.copy(),np.argwhere(isRef), useRatio=useRatio) + avNeg = funcForAveraging(data[isRef],axis=0) else: diff = data + avNeg = np.zeros_like(avData) # normalize signal for laser intensity if provided if lpower is not None: @@ -116,10 +118,10 @@ def averageScanPoints(scan,data,isRef=None,lpower=None,useRatio=False,\ # store error of mean err[i] = noise/np.sqrt(shot_idx.sum()) - - ret = dict(scan=scan_pos,data=ret,dataUnmasked=ret.copy(),err=err,errUnmasked=err.copy(), - chi2_0=chi2_0,diffsInScanPoint=diffsInScanPoint,dataAbsAvAll=avData,dataAbsAvScanPoint=data_abs, - dataAbs=data.copy()) + ret = dict(scan=scan_pos,data=ret,dataUnmasked=ret.copy(),err=err, + errUnmasked=err.copy(),chi2_0=chi2_0,diffsInScanPoint=diffsInScanPoint, + dataAbsAvNeg = avNeg, dataAsAbs=ret+avNeg, + dataAbsAvAll=avData,dataAbsAvScanPoint=data_abs,dataAbs=data.copy()) ret = storage.DataStorage(ret) return ret diff --git a/xray/dspace.py b/xray/dspace.py new file mode 100644 index 0000000..40591fd --- /dev/null +++ b/xray/dspace.py @@ -0,0 +1,70 @@ +from __future__ import division,print_function +import numpy as np +from numpy import sin,cos + + + +class Triclinic(object): + def __init__(self,a=1,b=1,c=1,alpha=90,beta=90,gamma=90): + self.a = a + self.b = b + self.c = c + alpha = alpha*np.pi/180 + beta = beta*np.pi/180 + gamma = gamma*np.pi/180 + self.alpha = alpha + self.beta = beta + self.gamma = gamma + + self._s11 = b**2 * c**2 * sin(alpha)**2 + self._s22 = a**2 * c**2 * sin(beta)**2 + self._s33 = a**2 * b**2 * sin(gamma)**2 + self._s12 = a*b*c**2*(cos(alpha) * cos(beta) - cos(gamma)) + self._s23 = a**2*b*c*(cos(beta) * cos(gamma) - cos(alpha)) + self._s13 = a*b**2*c*(cos(gamma) * cos(alpha) - cos(beta)) + self.V = (a*b*c)*np.sqrt(1-cos(alpha)**2 - cos(beta)**2 - cos(gamma)**2 + 2*cos(alpha)*cos(beta)*cos(gamma)) + + def __call__(self,h,k,l): return self.q(h,k,l) + + def d(self,h,k,l): + temp = self._s11*h**2 + \ + self._s22*k**2 + \ + self._s33*l**2 + \ + 2*self._s12*h*k+ \ + 2*self._s23*k*l+ \ + 2*self._s13*h*l + d = self.V/np.sqrt(temp) + return d + + def q(self,h,k,l): + return 2*np.pi/self.d(h,k,l) + +class Orthorombic(Triclinic): + def __init__(self,a=1,b=1,c=1): + Triclinic.__init__(self,a=a,b=b,c=c,alpha=90,beta=90,gamma=90) + +class Monoclinic(object): + def __init__(self,a=1,b=1,c=1,beta=90.): + self.a = a + self.b = b + self.c = c + beta = beta/np.pi*180 + self.beta = beta + + self.V = (a*b*c) + + def __call__(self,h,k,l): return self.Q(h,k,l) + + def Q(self,h,k,l): + temp = h**2/self.a**2 + (k*sin(self.beta))**2/self.b**2+l**2/self.c**2+2*h*l*cos(self.beta)/self.a/self.c + d = 1/np.sqrt(temp) + print(d) + return 2*np.pi/d + + + +ti3o5_lambda = Triclinic(a = 9.83776, b = 3.78674, c = 9.97069, beta = 91.2567) +ti3o5_beta = Triclinic(a = 9.7382 , b = 3.8005 , c = 9.4333 , beta = 91.496) +#ti3o5_beta = Monoclinic(a = 9.7382 , b = 3.8005 , c = 9.4333 , beta = 91.496) +ti3o5_alpha = Triclinic(a = 9.8372, b = 3.7921, c = 9.9717) +#ti3o5_alpha1 = Orthorombic(a = 9.8372, b = 3.7921, c = 9.9717) diff --git a/xray/id9.py b/xray/id9.py index cf29abf..5818f68 100644 --- a/xray/id9.py +++ b/xray/id9.py @@ -45,18 +45,21 @@ def doFolder_azav(folder,nQ=1500,force=False,mask=None,saveChi=True, return azav.doFolder(folder,files="*.edf*",nQ=nQ,force=force,mask=mask, saveChi=saveChi,poni=poni,storageFile=storageFile,diagnostic=diag) -def doFolder_dataRed(folder,storageFile='auto',monitor=None, - funcForAveraging=np.nanmean,errMask=5,chi2Mask=2,qlims=None): - """ storageFile cab the the basename of the file (to look for in folder) or - a DataStorage instance """ +def doFolder_dataRed(azavStorage,monitor=None,funcForAveraging=np.nanmean, + errMask=5,chi2Mask=2,qlims=None,outStorageFile='auto'): + """ azavStorage if a DataStorage instance or the filename to read """ - if storageFile == 'auto' : storageFile = folder + "/" + "pyfai_1d" + default_extension - - if isinstance(storageFile,storage.DataStorage): - data = storageFile + if isinstance(azavStorage,storage.DataStorage): + data = azavStorage + folder = os.path.dirname(data.filename) if data.filename is not None else "./" + elif os.path.exists(azavStorage): + folder = os.path.dirname(azavStorage) + data = storage.DataStorage(azavStorage) else: - # read azimuthal averaged curves - data = storage.DataStorage(storageFile) + # assume is just a folder name + folder = azavStorage + azavStorage = folder + "/pyfai_1d" + default_extension + data = storage.DataStorage(azavStorage) if qlims is not None: idx = (data.q>qlims[0]) & (data.qnp.argmax(y)][0] + c = int( (p1+p2)/2 ) + dp = int( np.abs(p1-p2) ) + idx = slice(c-dp,c+dp) + x = x[idx] + y = y[idx] + pars = pv.guess(y,x=x) + ret = pv.fit(y,x=x,params=pars) + return ret diff --git a/xray/storage.py b/xray/storage.py index 91430e8..5d41e2f 100644 --- a/xray/storage.py +++ b/xray/storage.py @@ -27,7 +27,7 @@ def unwrapArray(a,recursive=True,readH5pyDataset=True): a = dict(a); # convert to dict, otherwise can't asssign values for key,value in a.items(): a[key] = unwrapArray(value) elif isinstance(a,list): - for index in range(len(a)): a[index] = unwrapArray(a[i]) + for index in range(len(a)): a[index] = unwrapArray(a[index]) else: pass if isinstance(a,dict): a = DataStorage(a) @@ -80,25 +80,29 @@ def read(fname): extension = os.path.splitext(fname)[1] log.info("Reading storage file %s"%fname) if extension == ".npz": - return npzToDict(fname) + return DataStorage(npzToDict(fname)) elif extension == ".h5": - return h5ToDict(fname) + return DataStorage(h5ToDict(fname)) else: raise ValueError("Extension must be h5 or npz, it was %s"%extension) def save(fname,d): extension = os.path.splitext(fname)[1] log.info("Saving storage file %s"%fname) - if extension == ".npz": - return dictToNpz(fname,d) - elif extension == ".h5": - return dictToH5(fname,d) - else: - raise ValueError("Extension must be h5 or npz") + try: + if extension == ".npz": + return dictToNpz(fname,d) + elif extension == ".h5": + return dictToH5(fname,d) + else: + raise ValueError("Extension must be h5 or npz") + except Exception as e: + log.exception("Could not save %s"%fname) class DataStorage(dict): """ Storage for 1d integrated info """ - def __init__(self,fileOrDict,default_name='pyfai_1d',default_ext='npz'): + def __init__(self,fileOrDict,recursive=True, + default_name='pyfai_1d',default_ext='npz'): if isinstance(fileOrDict,dict): self.filename = None d = fileOrDict @@ -109,6 +113,11 @@ class DataStorage(dict): self.filename = fileOrDict d = read(fileOrDict) + if recursive: + for k in d.keys(): + if not isinstance(d[k],DataStorage) and isinstance(d[k],dict): + d[k] = DataStorage(d[k]) + # allow accessing with .data, .delays, etc. for k,v in d.items(): setattr(self,k,v) @@ -128,6 +137,31 @@ class DataStorage(dict): delattr(self,key) super().__delitem__(key) + def __str__(self): + keys = list(self.keys()) + keys.sort() + return "DataStorage obj containing: %s" % ",".join(keys) + + def __repr__(self): + keys = list(self.keys()) + keys.sort() + nchars = max(map(len,keys)) + fmt = "%%%ds %%s" % (nchars) + s = ["DataStorage obj containing (sorted): ",] + for k in keys: + if isinstance(self[k],np.ndarray): + value_str = "array %s"% "x".join(map(str,self[k].shape)) + elif isinstance(self[k],DataStorage): + value_str = str(self[k])[:50] + "..." + elif isinstance(self[k],(str,DataStorage)): + value_str = self[k][:50] + "..." + elif self[k] is None: + value_str = "None" + else: + value_str = str(self[k]) + s.append( fmt % (k,value_str) ) + return "\n".join(s) + def save(self,fname=None): if fname is None: fname = self.filename assert fname is not None diff --git a/xray/utils.py b/xray/utils.py index 552c409..4053f5c 100644 --- a/xray/utils.py +++ b/xray/utils.py @@ -9,6 +9,7 @@ import os import glob import pathlib import re +import numbers from . import storage as storage try: @@ -72,6 +73,30 @@ def removeExt(fname): def getBasename(fname): return removeExt(os.path.basename(fname)); +def findSlice(array,lims): + start = np.ravel(np.argwhere(array>lims[0]))[0] + stop = np.ravel(np.argwhere(array