From 93a1501a7be7fc85db07374a905465f2c4f4cac2 Mon Sep 17 00:00:00 2001 From: Thierry Zamofing Date: Tue, 8 Nov 2022 07:21:52 +0100 Subject: [PATCH] wip --- geometry.py | 204 ++++++++++++------------------------------------ graphExample.py | 107 ++++++++++++++++++------- 2 files changed, 129 insertions(+), 182 deletions(-) diff --git a/geometry.py b/geometry.py index b1435f4..e2e7f98 100644 --- a/geometry.py +++ b/geometry.py @@ -24,43 +24,30 @@ from scipy.optimize import fsolve class VLSgrating: def __init__(self): - # NEED TO CALL SETUP for parameters: - #self._param={'R':65.144m , # VLS curcature radius - # 'a0':1160.759, # gr/mm (gratings per mm) - # 'a1':0.3409, # gr/mm - # 'a2':-3.74E-5., # gr/mm^2 - # 'a3':1.98E-7 # gr/mm^3 - # } #spec - # gratings/mm = N(y)=a0+a1*y+a2*y^2+a3*y^3 - # - #self._lut_energy2geometry - #lookuptable for starting values. The lut is a numpy array with collums: - #energy (sort asc), R1, R2, alpha, beta, gamma' + # NEED TO CALL SETUP for parameters: _param and _equ pass - def setup(R,a0,a1,a2,a3,lut): - self._param={'R':R,'a0':a0,'a1':a1,'a2':a2,'a3':a3} - self._lut_energy2geometry=np.array(lut) - pass + def setup(self,vlsType='80meV'): + if vlsType=='80meV': + self._param={'R':65144.054, # radius of the grating mirror + 'a0':1160.759,'a1':0.3409,'a2':-3.74E-5,'a3':1.98E-7} # polynome for line density of the mirror + self._equ={'k': 1, # deflection order (default=1) + 'lut': # lookup table for starting guess values at a given energy for equation solver + [ # energy (eV), r1(mm), r2(mm), alpha(rad), + (200, 1650, 3500, np.deg2rad(88)), ] } + else: + raise ValueError('unknown type') - def setup_npz(self, fn): - #setup using a numpy npz file - pass - - def save_npz(self, fn): - #save current parameters as an numpy npz file - pass - - def update_lut(self, energyLst): - # recalculates the energy2geometry lut - # with the given energies - pass - - def energy2geometry(self, energy,mode,prec): + def energy2geometry(self, energy): # mode raytrace or interpolate lut # prec precision requitements (precision iteration) - # returns R1,R2,aa,gg - pass + # returns R1,R2,aa,bb,cc + pn=self.solve_focus_equ(energy) + r1, r2, alpha=pn + beta=self.beta(*pn) + gamma=self.gamma(*pn) + geo=(r1,r2,np.rad2deg(alpha), np.rad2deg(beta), np.rad2deg(gamma)) + return geo def geometry2raw(self, meas, debug=False): # returns raw motor positions @@ -68,126 +55,54 @@ class VLSgrating: # ---------- eugenio calculation functions ---------- #def bring_to_focus(self, eV): - def solve_focus_equ(self, eV,p0=None): + def solve_focus_equ(self, eV): es=self._equ - L, r1, alpha=es['L'], es['r1'], es['alpha'] # Take the current values from a dictionary - r2=L-r1 # starting R2 is taken as the difference between R1 and the maximum allowed length L - # p0=(np.deg2rad(88), 1650, 3500) - if not p0: - p0=[alpha, r1, r2] - p0=(np.deg2rad(88), 1650, 3500) - - + p0=es['lut'][0][1:] # fsolve starting values for [r1, r2, alpha] es['En']=eV # Update the dictionary with the energy at which I want to calculate the positions es['Lam_nm']=1239.842/eV # Convert energy to wavelength pn=fsolve(self.equations, p0, maxfev=100000000) - es['alpha'], es['r1'], es['r2']=pn + #es['r1'], es['r2'], es['alpha']=pn return pn def equations(self, pn): # system of 3 equations for 3 parameters return (self.Eq1(*pn), self.Eq2(*pn), self.Eq4(*pn)) - def Eq1(self, alpha, r1, r2): + def Eq1(self, r1, r2, alpha): p=self._param;equ=self._equ a1, R, a0, k, Lam_nm=p['a1'], p['R'], p['a0'], equ['k'], equ['Lam_nm'] - return np.cos(alpha)**2/r1+np.cos(self.beta(alpha, a0))**2/r2-\ - (np.cos(alpha)+np.cos(self.beta(alpha, a0)))/R-a1*k*Lam_nm/1e6 + beta=self.beta(r1, r2, alpha) + return np.cos(alpha)**2/r1+np.cos(beta)**2/r2-\ + (np.cos(alpha)+np.cos(beta))/R-a1*k*Lam_nm/1e6 - def Eq2(self, alpha, r1, r2): + def Eq2(self, r1, r2, alpha): p=self._param;equ=self._equ a0, a1, a2, R, k, Lam_nm=p['a0'], p['a1'], p['a2'], p['R'], equ['k'], equ['Lam_nm'] + beta=self.beta(r1, r2, alpha) return np.sin(alpha)/r1*(np.cos(alpha)**2/r1-np.cos(alpha)/R)-\ - np.sin(self.beta(alpha, a0))/r2*( - np.cos(self.beta(alpha, a0))**2/r2-np.cos(self.beta(alpha, a0))/R)+2*a2*k*Lam_nm/1e6/3 + np.sin(beta)/r2*( + np.cos(beta)**2/r2-np.cos(beta)/R)+2*a2*k*Lam_nm/1e6/3 - def Eq4(self, alpha, r1, r2): + def Eq4(self, r1, r2, alpha): p=self._param;equ=self._equ a3, R, k, Lam_nm, a0=p['a3'], p['R'], equ['k'], equ['Lam_nm'], p['a0'] + beta=self.beta(r1, r2, alpha) return 4*np.sin(alpha)**2/r1**2*(np.cos(alpha)**2/r1-np.cos(alpha)/R)\ -1/r1*(np.cos(alpha)**2/r1-np.cos(alpha)/R)**2+1/R**2*(1/r1-np.cos(alpha)/R)\ - +4*np.sin(self.beta(alpha, a0))**2/r2**2*( - np.cos(self.beta(alpha, a0))**2/r2-np.cos(self.beta(alpha, a0))/R)\ - -1/r2*(np.cos(self.beta(alpha, a0))**2/r2-np.cos(self.beta(alpha, a0))/R)**2\ - +1/R**2*(1/r2-np.cos(self.beta(alpha, a0))/R)-2*a3*k*Lam_nm/1e6 + +4*np.sin(beta)**2/r2**2*( + np.cos(beta)**2/r2-np.cos(beta)/R)\ + -1/r2*(np.cos(beta)**2/r2-np.cos(beta)/R)**2\ + +1/R**2*(1/r2-np.cos(beta)/R)-2*a3*k*Lam_nm/1e6 - def beta(self, alpha, a0): # calculate beta as a function of alpha, energy, and diffraction order (k) + def beta(self, r1, r2, alpha): # calculate beta as a function of alpha, energy, and diffraction order (k) p=self._param;equ=self._equ - k, Lam_nm=equ['k'], equ['Lam_nm'] + a0, k, Lam_nm=p['a0'], equ['k'], equ['Lam_nm'] return np.arcsin(np.sin(alpha)-a0*k*Lam_nm/1e6) - def chi(self, alpha, r1, r2): # compute detector lean angle + def gamma(self, r1, r2, alpha): # compute detector lean angle in respect of the r2 arm p=self._param;equ=self._equ - alpha, a0, r2, R, a1=equ['alpha'], p['a0'], equ['r2'], p['R'], p['a1'] - return np.arctan(np.cos(self.beta(alpha, a0))/(2*np.sin(self.beta(alpha, a0))-r2*(np.tan(self.beta(alpha, a0))/R+a1/a0))) - - -def testFunc(): - - def bring_to_focus(e, pars, std_guess): - L, r1, alpha=pars['L'], pars['r1'], pars['alpha'] # Take the current values from a dictionary - r2=L-r1 # starting R2 is taken as the difference between R1 and the maximum allowed length L - - if (abs(pars['En']-e)>200): # If the new energy is far from the current energy, use the user-defined starting guess - p0=std_guess - else: # otherwise, take the current values as starting guess - p0=[alpha, r1, r2] - - pars['En']=e # Update the dictionary with the energy at which I want to calculate the positions - pars['Lam_nm']=1239.842/e # Convert energy to wavelength - alpha, r1, r2=fsolve(equations, p0, maxfev=100000000, args=pars) - pars['alpha'], pars['r1'], pars['r2']=alpha, r1, r2 - return pars - - - def equations(p, pars): # system of 3 equations for 3 parameters - alpha, r1, r2=p - return (Eq1(alpha, r1, r2, pars), Eq2(alpha, r1, r2, pars), Eq4(alpha, r1, r2, pars)) - - - def Eq1(alpha, r1, r2, pars): - a1, R, a0, k, Lam_nm=pars['a1'], pars['R'], pars['a0'], pars['k'], pars['Lam_nm'] - return np.cos(alpha)**2/r1+np.cos(beta(alpha, a0, pars))**2/r2-\ - (np.cos(alpha)+np.cos(beta(alpha, a0, pars)))/R-a1*k*Lam_nm/1e6 - - - def Eq2(alpha, r1, r2, pars): - a0, a1, a2, R, k, Lam_nm=pars['a0'], pars['a1'], pars['a2'], pars['R'], pars['k'], pars['Lam_nm'] - return np.sin(alpha)/r1*(np.cos(alpha)**2/r1-np.cos(alpha)/R)-\ - np.sin(beta(alpha, a0, pars))/r2*( - np.cos(beta(alpha, a0, pars))**2/r2-np.cos(beta(alpha, a0, pars))/R)+2*a2*k*Lam_nm/1e6/3 - - - def Eq4(alpha, r1, r2, pars): - a3, R, k, Lam_nm, a0=pars['a3'], pars['R'], pars['k'], pars['Lam_nm'], pars['a0'] - return 4*np.sin(alpha)**2/r1**2*(np.cos(alpha)**2/r1-np.cos(alpha)/R)\ - -1/r1*(np.cos(alpha)**2/r1-np.cos(alpha)/R)**2+1/R**2*(1/r1-np.cos(alpha)/R)\ - +4*np.sin(beta(alpha, a0, pars))**2/r2**2*(np.cos(beta(alpha, a0, pars))**2/r2-np.cos(beta(alpha, a0, pars))/R)\ - -1/r2*(np.cos(beta(alpha, a0, pars))**2/r2-np.cos(beta(alpha, a0, pars))/R)**2\ - +1/R**2*(1/r2-np.cos(beta(alpha, a0, pars))/R)-2*a3*k*Lam_nm/1e6 - - - def beta(alpha, a0, pars): # calculate beta as a function of alpha, energy, and diffraction order (k) - k, Lam_nm=pars['k'], pars['Lam_nm'] - return np.arcsin(np.sin(alpha)-a0*k*Lam_nm/1e6) - - - def get_chi(pars): # compute detector lean angle - alpha, a0, r2, R, a1=pars['alpha'], pars['a0'], pars['r2'], pars['R'], pars['a1'] - return np.arctan( - np.cos(beta(alpha, a0, pars))/(2*np.sin(beta(alpha, a0, pars))-r2*(np.tan(beta(alpha, a0, pars))/R+a1/a0))) - - gratings=np.load('../doc/ForThierry/Standard.obj', allow_pickle='TRUE').item() - pars=gratings['Typ_80meV']['grating'] - - #Es = np.arange(250,1550, 50) - Es=[300,] - mypars = np.zeros(5) - for i,e in enumerate(Es): - pars = bring_to_focus(e,pars,[np.deg2rad(88),1650,3500]) - #pars = bring_to_focus(e,pars) - alpha, r1, r2, chi = pars['alpha'],pars['r1'],pars['r2'], get_chi(pars) - mypars = np.vstack((mypars, [e,np.rad2deg(alpha),r1,r2,np.rad2deg(chi)])) - print(mypars) + R, a0, a1=p['R'], p['a0'], p['a1'] + beta=self.beta(r1, r2, alpha) + return np.arctan(np.cos(beta)/(2*np.sin(beta)-r2*(np.tan(beta)/R+a1/a0))) if __name__=="__main__": import argparse @@ -229,52 +144,35 @@ if __name__=="__main__": np.set_printoptions(linewidth=196) if args.mode&0x01: - testFunc() gratings=np.load('../doc/ForThierry/Standard.obj', allow_pickle='TRUE').item() pars=gratings['Typ_80meV']['grating'] ### - - #self._param={'R':65.144,'a0':1160.759,'a1':0.3409,'a2':-3.74E-5.,'a3':1.98E-7} #spec - #self._param={'R':65.144,'a0':1160.77 ,'a1':0.3409,'a2':-3.73E-5.,'a3':1.98E-7} #production param - param={'R':65.144,'a0':1160.759,'a1':0.3409,'a2':-3.74E-5,'a3':1.98E-7, - 'lut': [#energy (eV), R1(mm), R2(mm), alpha(deg), beta(deg), gamma(deg) - [200 , 100., 300., .10, 20.0, 5.0], - [300 , 101., 301., .11, 20.1, 5.1], - [400 , 102., 302., .12, 20.2, 5.2], - [500 , 103., 303., .13, 20.3, 5.3], - ], - } param=dict() for k in ('R','a0','a1','a2','a3',): #fixed VLS parameters param[k]=pars.pop(k) #=pars[k] - - equ=dict() #equation solve parameters for k in ('L', 'r1', 'r2', 'alpha', 'k',): # starting guess values equ[k]=pars.pop(k) #=pars[k] - vlsg=VLSgrating() - vlsg._param=param - vlsg._equ=equ + vlsg.setup('80meV') + #vlsg._param=param + #vlsg._equ=equ + print(vlsg.energy2geometry(300)) Es=np.arange(250, 1550, 50) - Es=[300, ] + #Es=[300, ] mypars=[] for i, eV in enumerate(Es): pn=vlsg.solve_focus_equ(eV) - chi=vlsg.chi(*pn) - alpha, r1, r2=pn - mypars.append((eV, np.rad2deg(alpha), r1, r2, np.rad2deg(chi))) + r1, r2, alpha=pn + beta=vlsg.beta(*pn) + gamma=vlsg.gamma(*pn) + mypars.append((eV, r1, r2, np.rad2deg(alpha), np.rad2deg(beta), np.rad2deg(gamma))) mypars=np.array(mypars) print(mypars) - for eV in range(300,1500,50): - alpha, r1, r2 = vlsg.solve_focus_equ(eV) - - - #vlsg.setup(param) print(json.dumps(param, cls=MyJsonEncoder)) diff --git a/graphExample.py b/graphExample.py index 5b2089e..72938f1 100755 --- a/graphExample.py +++ b/graphExample.py @@ -16,6 +16,7 @@ import PyQt5.QtCore as QtCore import PyQt5.QtWidgets as QtW from PyQt5.uic import loadUiType import numpy as np +import geometry import sys @@ -130,12 +131,14 @@ class RIXSgrating(QWidget): super().__init__() self._param={'ctr':(400,400), # location of vlsg - 'aa':22, #grating angle - 'cc':58, #detector angle 'r1':314,#distance probe grating 'r2':447,#distance grating detector + 'aa':22, #grating angle + 'cc':58, #detector angle + #'geo':(8,8,5,5,1), # scaling: r1,r2,aa,bb,cc 'szG':(200,5), #size VLS grating 'szD':(150,5), #size detector + 'energy':300, # input energy in eV } self.initUI() @@ -152,26 +155,41 @@ class RIXSgrating(QWidget): self._wdGrpEnergy=w=QtGui.QGroupBox("Energy",self) w.move(10,10) #w.setFixedSize(200,42+32*2) - l=QtGui.QGridLayout(w) + lv=QtGui.QVBoxLayout(w) + w=QtGui.QWidget() #QGroupBox("E2",self) + lv.addWidget(w) + lg=QtGui.QGridLayout(w) + for i,t in enumerate(('Energy',)): w=QLabel(t) - l.addWidget(w, i,0) - w=QtGui.QLineEdit(t) - l.addWidget(w, i,1) + lg.addWidget(w, i,0) + w=QtGui.QLineEdit(t,objectName=t) + lg.addWidget(w, i,1) - self._wdGrpEnergy w=QtGui.QPushButton('calc geometry') - l.addWidget(w, i+1, 1) + i+=1;lg.addWidget(w, i, 1) - w=QtGui.QGroupBox("Geometry",self) - w.move(240,10) + sld={} + for key,rng,tk,pos in (('energy',(200,1500),50,60),): + sl=QtGui.QSlider(QtCore.Qt.Horizontal) + sl.setFixedWidth(200);sl.setMinimum(rng[0]);sl.setMaximum(rng[1]) + sl.setValue(self._param[key]) + sl.setTickPosition(QtGui.QSlider.TicksBelow);sl.setTickInterval(tk) + lv.addWidget(sl) + sl.valueChanged.connect(lambda val,key=key: self.sldChanged(key,val)) + sld[key]=sl + + + self._wdGrpGeometry=w=QtGui.QGroupBox("Geometry",self) + w.move(250,10) #w.setFixedSize(200,42+32*4) l=QtGui.QGridLayout(w) - #for i,t in enumerate(('R1','R2','alpha','beta','gamma',)): - for i, t in enumerate(('R1', 'R2', '\u03B1','\u03B2', '\u03B3',)): - w=QLabel(t) + lut={'aa':'\u03B1', 'bb':'\u03B2', 'cc':'\u03B3'} + for i,t in enumerate(('R1','R2','aa','bb','cc',)): + tl=lut.get(t,t) + w=QLabel(tl) l.addWidget(w, i,0) - w=QtGui.QLineEdit(t) + w=QtGui.QLineEdit(tl,objectName=t) l.addWidget(w, i,1) w=QtGui.QPushButton('calc raw motors') l.addWidget(w, i+1, 1) @@ -179,7 +197,7 @@ class RIXSgrating(QWidget): w=QtGui.QGroupBox("Drawing",self) - w.move(10,150) + w.move(10,160) l=QtGui.QVBoxLayout(w) sld={} @@ -196,9 +214,37 @@ class RIXSgrating(QWidget): def sldChanged(self,key,val,*args,**kwargs): print(key,val) - self._param[key]=val - self.update() + if key=='energy': + app=QApplication.instance() + wEnergy=self._wdGrpEnergy.findChild(QtGui.QLineEdit,'Energy') + wGeo=self._wdGrpGeometry + wR1=wGeo.findChild(QtGui.QLineEdit,'R1') + wR2=wGeo.findChild(QtGui.QLineEdit,'R2') + wAA=wGeo.findChild(QtGui.QLineEdit,'aa') + wBB=wGeo.findChild(QtGui.QLineEdit,'bb') + wCC=wGeo.findChild(QtGui.QLineEdit,'cc') + wEnergy.setText(f'{val:.4g}') + vlsg=app._vlsg + (r1,r2,aa,bb,cc)=vlsg.energy2geometry(val) + wR1.setText(f'{r1:.4g}') + wR2.setText(f'{r2:.4g}') + wAA.setText(f'{aa:.4g}') + wBB.setText(f'{bb:.4g}') + wCC.setText(f'{cc:.4g}') + p=self._param + + p['r1']=r1 + p['r2']=r2 + p['aa']=(90-aa) + p['bb']=(90-bb) + p['cc']=(90-cc) + else: + self._param[key]=val + + #print(v) + + self.update() def paintEvent(self, e): p=self._param @@ -211,6 +257,8 @@ class RIXSgrating(QWidget): szD=p['szD'] ctr=p['ctr'] + #scl=p['scl'] + qp = QPainter() qp.begin(self) @@ -274,7 +322,7 @@ class RIXSgrating(QWidget): col=QtGui.QColor() for i in range(n): p2=QtCore.QPointF(*tf3.map(szD[0]/2-szD[0]*i/(n-1),0)) - col.setHsv(i*scl,255,255,alpha) + col.setHsv(int(i*scl),255,255,alpha) pen.setColor(col) qp.setPen(pen) qp.drawLine(p0,p2) @@ -319,13 +367,13 @@ class RIXSgrating(QWidget): p9=QtCore.QPointF(*tf1.map(0, -w)) qp.drawLine(p2,p9) s=50 - qp.drawArc(int(p2.x())-s, int(p2.y())-s, 2*s, 2*s, 180*16, (aa)*16) + qp.drawArc(int(p2.x())-s, int(p2.y())-s, 2*s, 2*s, 180*16, int(aa*16)) s1=s*np.cos(aa*np.pi/180); s2=s*np.sin(aa*np.pi/180) qp.drawText(p2+QtCore.QPointF(-s1,s2+20),f'aa={aa:.5g}°') - qp.drawArc(int(p2.x())-s, int(p2.y())-s, 2*s, 2*s, aa*16, bb*16) + qp.drawArc(int(p2.x())-s, int(p2.y())-s, 2*s, 2*s, int(aa*16), int(bb*16)) qp.drawText(p2+QtCore.QPointF(s1,-s2+20),f'bb={bb:.5g}°') - qp.drawArc(int(p6.x())-s, int(p6.y())-s, 2*s, 2*s, (180+aa+bb)*16, (cc)*16) + qp.drawArc(int(p6.x())-s, int(p6.y())-s, 2*s, 2*s, int(180+aa+bb*16), int(cc*16)) p10=QtCore.QPointF(*tf2.map(r2+20, 0)) qp.drawText(p10,f'cc={cc:.5g}°') @@ -336,13 +384,14 @@ class RIXSgrating(QWidget): - -def main(): - app = QApplication(sys.argv) - #ex = RIXSgirder() - ex = RIXSgrating() +if __name__ == '__main__': + def main(): + app=QApplication(sys.argv) + vlsg=geometry.VLSgrating() + vlsg.setup('80meV') + app._vlsg=vlsg + # ex = RIXSgirder() + ex=RIXSgrating() sys.exit(app.exec_()) - -if __name__ == '__main__': - main() \ No newline at end of file + main() \ No newline at end of file