Furka/EsfRixsApps
0
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

wip

Browse Source
This commit is contained in:
zamofing_t
2022-11-08 07:21:52 +01:00
parent 6b9c6e5611
commit 93a1501a7b
2 changed files with 129 additions and 182 deletions
Download Patch File Download Diff File Expand all files Collapse all files

204
geometry.py
View File

@@ -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))

107
graphExample.py
View File

@@ -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()
main()