various changes: some naming, added function for averaging and changed doShots to accept list of calibs
This commit is contained in:
Marco Cammarata
parent
b0914f15f8
commit
a2c508fa5d
|
|
@ -101,6 +101,54 @@ def showShots(im1,im2):
|
|||
a[0].plot(p1)
|
||||
a[1].plot(p2)
|
||||
|
||||
def ratioOfAverage(p1,p2,threshold=0.03):
|
||||
"""
|
||||
p1 and p2 are the energy spectrum. if 2D the first index has to be the shot number
|
||||
calculate median ratio taking into account only regions where p1 and p2 are > 5% of the max """
|
||||
# check if they are 2D
|
||||
if p1.ndim == 1:
|
||||
p1 = p1[np.newaxis,:]
|
||||
p2 = p2[np.newaxis,:]
|
||||
# w1 and w2 are the weights
|
||||
w1 = p1.copy(); w2 = p2.copy()
|
||||
if threshold is not None:
|
||||
# weights will be set to zero if intensity is smaller than 5% of max
|
||||
# for each shots, get maximum
|
||||
m1 = np.nanmax(p1,axis=1); m2 = np.nanmax(p2,axis=1)
|
||||
# find where each spectrum is smaller than threshold*max_for_that_shot; they will be masked out
|
||||
idx1 = p1 < (m1[:,np.newaxis]*threshold)
|
||||
idx2 = p2 < (m2[:,np.newaxis]*threshold)
|
||||
w1[idx1]=0
|
||||
w2[idx2]=0
|
||||
# using masked array because some pixel will have zero shots contributing
|
||||
av1 = np.ma.average(p1,axis=0,weights=w1)
|
||||
av1[av1.mask] = np.nan
|
||||
av2 = np.ma.average(p2,axis=0,weights=w2)
|
||||
av2[av2.mask] = np.nan
|
||||
return av2/av1
|
||||
|
||||
def medianRatio(p1,p2,threshold=0.03):
|
||||
"""
|
||||
p1 and p2 are the energy spectrum. if 2D the first index has to be the shot number
|
||||
calculate median ratio taking into account only regions where p1 and p2 are > 5% of the max """
|
||||
# check if they are 2D
|
||||
if p1.ndim == 1:
|
||||
p1 = p1[np.newaxis,:]
|
||||
p2 = p2[np.newaxis,:]
|
||||
p1 = np.ma.asarray( p1.copy() )
|
||||
p2 = np.ma.asarray( p2.copy() )
|
||||
if threshold is not None:
|
||||
m1 = np.nanmax(p1,axis=1); m2 = np.nanmax(p2,axis=1)
|
||||
# find where each spectrum is smaller than threshold*max_for_that_shot; they will be masked out
|
||||
idx1 = p1 < (m1[:,np.newaxis]*threshold)
|
||||
idx2 = p2 < (m2[:,np.newaxis]*threshold)
|
||||
idx = idx1 & idx2
|
||||
p1.mask = idx
|
||||
p2.mask = idx
|
||||
ratio = p2/p1
|
||||
return np.ma.average(ratio,axis=0,weights=p1)
|
||||
|
||||
|
||||
|
||||
class AnalyzeRun(object):
|
||||
def __init__(self,run,initAlign="auto",swapx=g_swapx,swapy=g_swapy):
|
||||
|
|
@ -126,6 +174,7 @@ class AnalyzeRun(object):
|
|||
d = self.data
|
||||
self.spec1 = d.spec1 ; # spec1 is the one that is moved
|
||||
self.spec2 = d.spec2 ;
|
||||
self.E = alignment.defaultE
|
||||
|
||||
try:
|
||||
self.loadTransform(initAlign)
|
||||
|
|
@ -136,10 +185,15 @@ class AnalyzeRun(object):
|
|||
else:
|
||||
self.initAlign = initAlign
|
||||
|
||||
def getShot(self,shot=0,calib=None,bkgSub="line",roi=g_roi_height):
|
||||
def getShots(self,shots=0,calib=None,bkgSub="line",roi=g_roi_height):
|
||||
if shots == "all":
|
||||
if calib != None:
|
||||
shots = slice(0,self.nShotsPerCalib[calib])
|
||||
else:
|
||||
shots = slice(0,self.data.spec1.nShots)
|
||||
# read data
|
||||
im1 = self.spec1.getShots(shot,calib=calib)
|
||||
im2 = self.spec2.getShots(shot,calib=calib)
|
||||
im1 = self.spec1.getShots(shots,calib=calib)
|
||||
im2 = self.spec2.getShots(shots,calib=calib)
|
||||
# subtractBkg bkg
|
||||
im1 = alignment.subtractBkg(im1,bkg_type=bkgSub)
|
||||
im2 = alignment.subtractBkg(im2,bkg_type=bkgSub)
|
||||
|
|
@ -174,31 +228,40 @@ class AnalyzeRun(object):
|
|||
gui.save(fname)
|
||||
|
||||
|
||||
def analyzeScan(self,initpars=None,nShotsPerCalib=20,nC=None,doFit=False,fitEveryCalib=False,nSaveImg=5):
|
||||
""" this is a comment """
|
||||
def analyzeScan(self,initpars=None,shots=slice(0,30),calibs="all",nImagesToFit=0,nSaveImg=5):
|
||||
""" nImagesToFit: number of images to Fit per calibcycle, (int or "all") """
|
||||
if initpars is None: initpars= self.initAlign
|
||||
if nC is None: nC = self.nCalib
|
||||
out = []
|
||||
for i in range(nC):
|
||||
if nShotsPerCalib == 'all':
|
||||
shots = slice(self.nShotsPerCalib[i])
|
||||
if calibs == "all": calibs=list(range(self.nCalib))
|
||||
if isinstance(calibs,slice): calibs=list(range(self.nCalib))[calibs]
|
||||
nC = len(calibs)
|
||||
for ic,calib in enumerate(calibs):
|
||||
if nImagesToFit == "all":
|
||||
nToFit = self.nShotsPerCalib[calib]
|
||||
else:
|
||||
shots = slice(nShotsPerCalib)
|
||||
s1,s2 = self.getShot(shots,calib=i)
|
||||
if fitEveryCalib is not False:
|
||||
ret,bestTransf = alignment.doShots(s1[:fitEveryCalib],s2[:fitEveryCalib],doFit=True,\
|
||||
nToFit = nImagesToFit
|
||||
#print("Memory available 1",x3py.toolsOS.memAvailable())
|
||||
s1,s2 = self.getShots(shots,calib=calib)
|
||||
#print("Memory available 2",x3py.toolsOS.memAvailable())
|
||||
if nImagesToFit > 0:
|
||||
ret,bestTransf = alignment.doShots(s1[:nToFit],s2[:nToFit],doFit=True,\
|
||||
initpars=initpars,nSaveImg=nSaveImg,returnBestTransform=True);
|
||||
initpars = bestTransf; self.initAlign=bestTransf
|
||||
ret = alignment.doShots(s1,s2,initpars=initpars,doFit=doFit,nSaveImg=nSaveImg)
|
||||
self.results[i] = ret
|
||||
print("Calib cycle %d -> %.3f (best FOM: %.2f)" % (i,self.scanpos[i],np.min(ret.fom)))
|
||||
out.append(ret)
|
||||
return out
|
||||
if nToFit < s1.shape[0]:
|
||||
ret2 = alignment.doShots(s1[nToFit:],s2[nToFit:],initpars=initpars,doFit=False,nSaveImg=0)
|
||||
if ret is None:
|
||||
ret = ret2
|
||||
else:
|
||||
ret = alignment.unravel_results( (ret,ret2) )
|
||||
#print("Memory available 3",x3py.toolsOS.memAvailable())
|
||||
self.results[calib] = ret
|
||||
#print("Memory available 4",x3py.toolsOS.memAvailable())
|
||||
print("Calib cycle %d/%d -> %.3f (best FOM: %.2f)" % (ic,nC,self.scanpos[ic],np.nanmin(ret.fom)))
|
||||
return [self.results[c] for c in calibs]
|
||||
|
||||
def doShot(self,shot=0,calib=None,initpars=None,im1=None,im2=None,doFit=True,show=False,showInit=False,save=False,savePlot="auto"):
|
||||
if initpars is None: initpars= self.initAlign
|
||||
if (im1 is None) or (im2 is None):
|
||||
im1,im2 = self.getShot(shot,calib=calib); im1=im1[0]; im2=im2[0]
|
||||
im1,im2 = self.getShots(shot,calib=calib); im1=im1[0]; im2=im2[0]
|
||||
r = alignment.doShot(im1,im2,initpars,doFit=doFit,show=showInit)
|
||||
im1 = r.im1
|
||||
im2 = r.im2
|
||||
|
|
@ -219,7 +282,7 @@ class AnalyzeRun(object):
|
|||
"""
|
||||
if initpars is None: initpars= self.initAlign
|
||||
if shots == "all": shots = slice(self.nShotsPerCalib[calib])
|
||||
s1,s2 = self.getShot(shots,calib=calib)
|
||||
s1,s2 = self.getShots(shots,calib=calib)
|
||||
ret,transformForBestFit = alignment.doShots(s1,s2,initpars=initpars,doFit=doFit,\
|
||||
returnBestTransform=True,nSaveImg=nSaveImg)
|
||||
if doFit: self.initAlign = transformForBestFit
|
||||
|
|
@ -234,48 +297,34 @@ class AnalyzeRun(object):
|
|||
if len(self.results) == 0: print("self.results are empty, returning without saving")
|
||||
if not os.path.isdir(g_folder_out): os.makedirs(g_folder_out)
|
||||
if fname == "auto":
|
||||
fname = g_folder_out+"/run%04d_analysis.h5" % self.run
|
||||
fname = g_folder_out+"/run%04d_analysis" % self.run
|
||||
if os.path.exists(fname) and not overwrite:
|
||||
print("File %s exists, **NOT** saving, use overwrite=True is you want ..."%fname)
|
||||
return
|
||||
if os.path.exists(fname) and overwrite: os.unlink(fname)
|
||||
print("Saving results to %s"%fname)
|
||||
h = h5py.File(fname)
|
||||
h["roi1"] = (self.roi1.start,self.roi1.stop)
|
||||
h["roi2"] = (self.roi2.start,self.roi2.stop)
|
||||
h["scanmot0"] = self.data.scan.scanmotor0
|
||||
h["scanpos0"] = self.data.scan.scanmotor0_values
|
||||
if hasattr(self.data.scan,"scanmotor1"):
|
||||
h["scanmot1"] = self.data.scan.scanmotor1
|
||||
h["scanpos1"] = self.data.scan.scanmotor1_values
|
||||
#h["transform"] = self.initAlign
|
||||
for (c,v) in self.results.items():
|
||||
cname = "calib%04d/" % c if isinstance(c,int) else "calib%s/" % c
|
||||
for p,vv in mc.objToDict(v).items():
|
||||
# cannot save in hfd5 certain python objects
|
||||
if p == "fit_result" or p.find("final_transform")==0:
|
||||
continue
|
||||
if isinstance(vv,dict):
|
||||
for pname,parray in vv.items():
|
||||
name = cname + p + "/" + pname
|
||||
h[name] = parray
|
||||
else:
|
||||
h[cname + p] = vv
|
||||
h.close()
|
||||
|
||||
def load(self,fname="auto"):
|
||||
if fname == "auto": fname = g_folder_out+"/run%04d_analysis.npz" % self.run
|
||||
temp = np.load(fname)
|
||||
self.results = temp["results"].item()
|
||||
temp.close()
|
||||
|
||||
def saveTransform(self,fname="auto",transform=None):
|
||||
def _auto_transform_name(self,run=None,calib=None):
|
||||
if run is None: run = self.run
|
||||
fname = g_folder_init+"/run%04d_transform" % run
|
||||
if calib is not None:
|
||||
fname = fname + "_c%03d" % calib
|
||||
return fname + ".npy"
|
||||
|
||||
def saveTransform(self,fname="auto",calib=None,transform=None):
|
||||
if transform is None: transform = self.initAlign
|
||||
if fname == "auto":
|
||||
fname = g_folder_init+"/run%04d_transform.npy" % self.run
|
||||
if fname == "auto": fname = self._auto_transform_name(calib=calib)
|
||||
print("Saving roi and transformation parameter to %s"%fname)
|
||||
alignment.saveAlignment(fname,self.initAlign,self.roi1,self.roi2)
|
||||
|
||||
def loadTransform(self,fname="auto"):
|
||||
def loadTransform(self,fname="auto", calib=None):
|
||||
if isinstance(fname,dict): raise FileNotFoundError
|
||||
if fname == "auto":
|
||||
fname = g_folder_init+"/run%04d_transform.npy" % self.run
|
||||
elif isinstance(fname,int):
|
||||
if fname == "auto": fname = self._auto_transform_name(calib=calib)
|
||||
if isinstance(fname,int):
|
||||
fname = g_folder_init+"/run%04d_transform.npy" % fname
|
||||
if not os.path.exists(fname): print("Asked to read %s, but it does not exist"%fname)
|
||||
temp = np.load(fname).item()
|
||||
|
|
|
|||
Loading…
Reference in New Issue