From 66ca11e0142f0899b497076dedf354f150652a04 Mon Sep 17 00:00:00 2001 From: Marco Cammarata Date: Thu, 1 Dec 2016 11:14:42 +0100 Subject: [PATCH] various improvements (i.e. i forgot ...) --- alignment.py | 37 +++++++++++++++--- xanes_analyzeRun.py | 91 +++++++++++++++++++++++---------------------- 2 files changed, 77 insertions(+), 51 deletions(-) diff --git a/alignment.py b/alignment.py index d1e6d3e..30f0bea 100644 --- a/alignment.py +++ b/alignment.py @@ -16,6 +16,7 @@ import time #defaultE = np.arange(1024)*0.12+7060 defaultE = (np.arange(1024)-512)*0.189+7123 +defaultE = 7063.5+0.1295*np.arange(1024); # done on nov 30 2016, based on xppl3716:r77 __i = np.arange(1024) __x = (__i-512)/512 @@ -83,7 +84,7 @@ fit_ret = collections.namedtuple("fit_ret",["init_pars","final_pars",\ def calcFOM(p1,p2,ratio): idx = ( p1>p1.max()/10 )# & (p2>p2.max()/10) ratio = ratio[idx] - return ratio.std()/ratio.mean() + return ratio.std()/np.abs(ratio.mean()) def subtractBkg(imgs,nPix=100,bkg_type="line"): if imgs.ndim == 2: imgs = imgs[np.newaxis,:] @@ -265,6 +266,26 @@ def saveAlignment(fname,transform,roi1,roi2): def loadAlignment(fname): return np.load(fname).item() +def getBestTransform(results): + # try to see if it is unravelled + if not hasattr(results,"_fields"): results = unravel_results(results) + # find best based on FOM + idx = np.nanargmin(np.abs(results.fom)) + parnames = results.init_pars.keys() + bestPars = dict() + for n in parnames: + bestPars[n] = results.final_pars[n][idx] + if n.find("limit_") == 0: + if bestPars[n][0] is None: + bestPars[n] = None + else: + bestPars[n] = tuple(bestPars[n]) + return bestPars,np.nanmin(np.abs(results.fom)) + + + + + def unravel_results(res,nSaveImg='all'): final_pars = dict() # res[0].fit_result is a list if we are trying to unravel retult from getShots @@ -275,8 +296,12 @@ def unravel_results(res,nSaveImg='all'): final_pars = dict() init_pars = dict() for n in parnames: - final_pars[n] = np.hstack( [r.final_pars[n] for r in res]) - init_pars[n] = np.hstack( [r.init_pars[n] for r in res]) + if n.find("limit_") == 0: + final_pars[n] = np.vstack( [r.final_pars[n] for r in res]) + init_pars[n] = np.vstack( [r.init_pars[n] for r in res] ) + else: + final_pars[n] = np.hstack( [r.final_pars[n] for r in res]) + init_pars[n] = np.hstack( [r.init_pars[n] for r in res] ) im1 = np.asarray( [r.im1 for r in res if r.im1.shape[0] != 0]) im2 = np.asarray( [r.im2 for r in res if r.im2.shape[0] != 0]) if nSaveImg != 'all': @@ -582,9 +607,9 @@ def doShots(imgs1,imgs2,initpars,nJobs=16,doFit=False,returnBestTransform=False, (joblib.delayed(doShot)(imgs1[i],imgs2[i],initpars,doFit=doFit) for i in range(N)) res = unravel_results(pool,nSaveImg=nSaveImg) if returnBestTransform: - idx = np.nanargmin( res.fom ); - print("FOM for best alignment %.2f"%res.fom[idx]) - return res,pool[idx].final_pars + bestT,fom = getBestTransform(res) + print("FOM for best alignment %.2f"%fom) + return res,bestT else: return res # out = collections.OrderedDict( enumerate(pool) ) diff --git a/xanes_analyzeRun.py b/xanes_analyzeRun.py index d055edc..7c61a64 100644 --- a/xanes_analyzeRun.py +++ b/xanes_analyzeRun.py @@ -221,6 +221,50 @@ class AnalyzeRun(object): self.initAlign = gui.transform gui.save(fname) + 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.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 + self.initAlign = r.final_pars + if show: + if savePlot == "auto": + if not os.path.isdir(g_folder_out): os.makedirs(g_folder_out) + savePlot = g_folder_out+"/run%04d_calib%s_shot%04d_fit.png" % (self.run,calib,shot) + alignment.plotShot(im1,im2,res=r,save=savePlot) + if save: self.saveTransform() + return r + + def doShots(self,shots=slice(0,50),calib=None,initpars=None,doFit=False,returnBestTransform=False,nSaveImg=5,nInChunks=250): + """ + shots : slice to define shots to read, use 'all' for all shots in calibcycle + nSaveImg : save saveImg images in memory (self.results), use 'all' for all + useful for decreasing memory footprint + """ + if initpars is None: initpars= self.initAlign + if shots == "all": shots = list(range(self.nShotsPerCalib[calib])) + if isinstance(shots,slice): + nmax = self.nShotsPerCalib[calib] if calib is not None else self.data.spec1.nShots + shots = sliceToIndices(shots,nmax) + chunks = x3py.toolsVarious.chunk(shots,nInChunks) + ret_chunks = [] + for chunk in chunks: + s1,s2 = self.getShots(chunk,calib=calib) + ret_chunk = alignment.doShots(s1,s2,initpars=initpars,doFit=doFit,\ + returnBestTransform=False,nSaveImg=nSaveImg) + ret_chunks.append(ret_chunk) + if nSaveImg != 'all' and len(chunk) %.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.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 - self.initAlign = r.final_pars - if show: - if savePlot == "auto": - if not os.path.isdir(g_folder_out): os.makedirs(g_folder_out) - savePlot = g_folder_out+"/run%04d_calib%s_shot%04d_fit.png" % (self.run,calib,shot) - alignment.plotShot(im1,im2,res=r,save=savePlot) - if save: self.saveTransform() - return r - - def doShots(self,shots=slice(0,50),calib=None,initpars=None,doFit=False,returnBestTransform=False,nSaveImg='all',nInChunks=250): - """ - shots : slice to define shots to read, use 'all' for all shots in calibcycle - nSaveImg : save saveImg images in memory (self.results), use 'all' for all - useful for decreasing memory footprint - """ - if initpars is None: initpars= self.initAlign - if shots == "all": shots = list(range(slice(self.nShotsPerCalib[calib]))) - if isinstance(shots,slice): - nmax = self.nShotsPerCalib[calib] if calib is not None else self.data.spec1.nShots - shots = sliceToIndices(shots,nmax) - chunks = x3py.toolsVarious.chunk(shots,nInChunks) - ret = [] - for chunk in chunks: - s1,s2 = self.getShots(chunk,calib=calib) - temp = alignment.doShots(s1,s2,initpars=initpars,doFit=doFit,\ - returnBestTransform=False,nSaveImg=nSaveImg) - ret.append(temp) - ret = alignment.unravel_results(ret) - idxBest = np.nanargmin(ret.fom) - transformForBestFit = dict( [ (p,ret.final_pars[p][idxBest]) for p in ret.final_pars.keys() ] ) - if doFit: self.initAlign = transformForBestFit - # keep it for later ! - self.results[calib] = ret - if returnBestTransform: - return ret,transformForBestFit - else: - return ret def save(self,fname="auto",overwrite=False): if len(self.results) == 0: