Update gaussianoptics.py
This commit is contained in:
parent
47ea699f3e
commit
b894e76e9f
@ -1,254 +1,254 @@
|
|||||||
import numpy as np
|
import numpy as np
|
||||||
|
|
||||||
def gaussianoptics(opticalbench,lam,q0,startposition,N):
|
def gaussianoptics(opticalbench,lam,q0,startposition,N):
|
||||||
'''
|
'''
|
||||||
'path', refractive index , position[mm], (relative: 'r', absolute: 'a' (position is the end position), filling: 'f')
|
'path', refractive index , position[mm], (relative: 'r', absolute: 'a' (position is the end position), filling: 'f')
|
||||||
'thinlens', f [mm] , position[mm], (relative: 'r', absolute: 'a')
|
'thinlens', f [mm] , position[mm], (relative: 'r', absolute: 'a')
|
||||||
'thicklens', [ n1 = refractive index outside of the lens,
|
'thicklens', [ n1 = refractive index outside of the lens,
|
||||||
n2 = refractive index of the lens itself (inside the lens),
|
n2 = refractive index of the lens itself (inside the lens),
|
||||||
R1 = Radius of curvature of First surface,
|
R1 = Radius of curvature of First surface,
|
||||||
R2 = Radius of curvature of Second surface,
|
R2 = Radius of curvature of Second surface,
|
||||||
t = center thickness of lens] , position[mm], (relative: 'r', absolute: 'a')
|
t = center thickness of lens] , position[mm], (relative: 'r', absolute: 'a')
|
||||||
'flatmirror', 0 , position[mm], (relative: 'r', absolute: 'a')
|
'flatmirror', 0 , position[mm], (relative: 'r', absolute: 'a')
|
||||||
'curvedmirror', [R [mm],
|
'curvedmirror', [R [mm],
|
||||||
angle [°],
|
angle [°],
|
||||||
plane (0=tangential, 1=sagittal)], position[mm], (relative: 'r', absolute: 'a')
|
plane (0=tangential, 1=sagittal)], position[mm], (relative: 'r', absolute: 'a')
|
||||||
'flatrefraction', [ n1 = initial refractive index,
|
'flatrefraction', [ n1 = initial refractive index,
|
||||||
n2 = final refractive index] , position[mm], (relative: 'r', absolute: 'a')
|
n2 = final refractive index] , position[mm], (relative: 'r', absolute: 'a')
|
||||||
'curvedrefraction', [ n1 = initial refractive index,
|
'curvedrefraction', [ n1 = initial refractive index,
|
||||||
n2 = final refractive index,
|
n2 = final refractive index,
|
||||||
R = radius of curvature R > 0 for convex (centre of
|
R = radius of curvature R > 0 for convex (centre of
|
||||||
curvature after interface)], position[mm], (relative: 'r', absolute: 'a')
|
curvature after interface)], position[mm], (relative: 'r', absolute: 'a')
|
||||||
example: opticalbench = ['path',1,0.050,'a''thinlens',0.05,0.100,'r''path',[1,0,0,0,0],0,'f''thinlens',[0.15,0,0,0,0],0.200,'r']
|
example: opticalbench = ['path',1,0.050,'a''thinlens',0.05,0.100,'r''path',[1,0,0,0,0],0,'f''thinlens',[0.15,0,0,0,0],0.200,'r']
|
||||||
TODO:
|
TODO:
|
||||||
- check if 'path' is first and last element
|
- check if 'path' is first and last element
|
||||||
- find 'empty space'
|
- find 'empty space'
|
||||||
- normalize curve by refractive index other than 1
|
- normalize curve by refractive index other than 1
|
||||||
Jan Arenskötter 05.08.2020
|
Elena Arenskötter 05.08.2020
|
||||||
|
|
||||||
'''
|
'''
|
||||||
|
|
||||||
elements = len(opticalbench)
|
elements = len(opticalbench)
|
||||||
OB = [[0,0,0,0] for j in range(elements)] #[optical element, property, start position, stop position]
|
OB = [[0,0,0,0] for j in range(elements)] #[optical element, property, start position, stop position]
|
||||||
position = startposition
|
position = startposition
|
||||||
|
|
||||||
#find absolute positions
|
#find absolute positions
|
||||||
for i in range(0,elements):
|
for i in range(0,elements):
|
||||||
if opticalbench[i][0] == 'path':
|
if opticalbench[i][0] == 'path':
|
||||||
OB[i][0] = 'path'
|
OB[i][0] = 'path'
|
||||||
OB[i][1] = opticalbench[i][1] #properties
|
OB[i][1] = opticalbench[i][1] #properties
|
||||||
#start position
|
#start position
|
||||||
OB[i][2] = position
|
OB[i][2] = position
|
||||||
#stop position
|
#stop position
|
||||||
if opticalbench[i][3] == 'r':
|
if opticalbench[i][3] == 'r':
|
||||||
position = position + opticalbench[i][2]
|
position = position + opticalbench[i][2]
|
||||||
OB[i][3] = position
|
OB[i][3] = position
|
||||||
elif opticalbench[i][3] == 'a':
|
elif opticalbench[i][3] == 'a':
|
||||||
if opticalbench[i][2] < position:
|
if opticalbench[i][2] < position:
|
||||||
disp('Check sequential arrangement!')
|
disp('Check sequential arrangement!')
|
||||||
break
|
break
|
||||||
else:
|
else:
|
||||||
position = opticalbench[i][2]
|
position = opticalbench[i][2]
|
||||||
OB[i][3] = position
|
OB[i][3] = position
|
||||||
|
|
||||||
elif opticalbench[i][3] == 'f':
|
elif opticalbench[i][3] == 'f':
|
||||||
if i == elements-1:
|
if i == elements-1:
|
||||||
position = position + 100
|
position = position + 100
|
||||||
OB[i][3] = position
|
OB[i][3] = position
|
||||||
print('Warning: Floating path at the end of your optical bench!')
|
print('Warning: Floating path at the end of your optical bench!')
|
||||||
print('-> Replaced by 100mm relative length.')
|
print('-> Replaced by 100mm relative length.')
|
||||||
print('Optical bench entry: '+str(i))
|
print('Optical bench entry: '+str(i))
|
||||||
|
|
||||||
else:
|
else:
|
||||||
if opticalbench[i+1][3] == 'a':
|
if opticalbench[i+1][3] == 'a':
|
||||||
position = opticalbench[i+1][2]
|
position = opticalbench[i+1][2]
|
||||||
OB[i][3] = position
|
OB[i][3] = position
|
||||||
elif opticalbench[i+1][3] == 'r':
|
elif opticalbench[i+1][3] == 'r':
|
||||||
position = position + opticalbench[i+1][2]
|
position = position + opticalbench[i+1][2]
|
||||||
OB[i][3] = position
|
OB[i][3] = position
|
||||||
elif opticalbench[i+1][3] == 'f':
|
elif opticalbench[i+1][3] == 'f':
|
||||||
print('Warning: Consecutive floating/relative elements!')
|
print('Warning: Consecutive floating/relative elements!')
|
||||||
print('Optical bench entry: '+str(i))
|
print('Optical bench entry: '+str(i))
|
||||||
break
|
break
|
||||||
else:
|
else:
|
||||||
print('Warning: Unknown position parameter.')
|
print('Warning: Unknown position parameter.')
|
||||||
print('Optical bench entry: ',str(i+1))
|
print('Optical bench entry: ',str(i+1))
|
||||||
break
|
break
|
||||||
|
|
||||||
else:
|
else:
|
||||||
print('Warning: Unknown position parameter.')
|
print('Warning: Unknown position parameter.')
|
||||||
print('Optical bench entry: '+str(i))
|
print('Optical bench entry: '+str(i))
|
||||||
break
|
break
|
||||||
|
|
||||||
else:
|
else:
|
||||||
if opticalbench[i][0] == 'thinlens':
|
if opticalbench[i][0] == 'thinlens':
|
||||||
OB[i][0] = 'thinlens'
|
OB[i][0] = 'thinlens'
|
||||||
elif opticalbench[i][0] == 'thicklens':
|
elif opticalbench[i][0] == 'thicklens':
|
||||||
OB[i][0] = 'thicklens'
|
OB[i][0] = 'thicklens'
|
||||||
elif opticalbench[i][0] == 'flatmirror':
|
elif opticalbench[i][0] == 'flatmirror':
|
||||||
OB[i][0] = 'flatmirror'
|
OB[i][0] = 'flatmirror'
|
||||||
elif opticalbench[i][0] == 'curvedmirror':
|
elif opticalbench[i][0] == 'curvedmirror':
|
||||||
OB[i][0] = 'curvedmirror'
|
OB[i][0] = 'curvedmirror'
|
||||||
elif opticalbench[i][0] == 'flatrefraction':
|
elif opticalbench[i][0] == 'flatrefraction':
|
||||||
OB[i][0] = 'flatrefraction'
|
OB[i][0] = 'flatrefraction'
|
||||||
elif opticalbench[i][0] == 'curvedrefraction':
|
elif opticalbench[i][0] == 'curvedrefraction':
|
||||||
OB[i][0] = 'curvedrefraction'
|
OB[i][0] = 'curvedrefraction'
|
||||||
else:
|
else:
|
||||||
print('Warning: Unknown optical element!')
|
print('Warning: Unknown optical element!')
|
||||||
print('--> '+opticalbench[i][0]+' entry: '+str(i))
|
print('--> '+opticalbench[i][0]+' entry: '+str(i))
|
||||||
break
|
break
|
||||||
|
|
||||||
OB[i][1]=opticalbench[i][1] #properties
|
OB[i][1]=opticalbench[i][1] #properties
|
||||||
#start position
|
#start position
|
||||||
OB[i][2] = position
|
OB[i][2] = position
|
||||||
#stop position
|
#stop position
|
||||||
if opticalbench[i][3] == 'a':
|
if opticalbench[i][3] == 'a':
|
||||||
if opticalbench[i][2] < position:
|
if opticalbench[i][2] < position:
|
||||||
print('Check sequential arrangement!')
|
print('Check sequential arrangement!')
|
||||||
print('Optical bench entry: '+str(i))
|
print('Optical bench entry: '+str(i))
|
||||||
break
|
break
|
||||||
else:
|
else:
|
||||||
position = opticalbench[i][2]
|
position = opticalbench[i][2]
|
||||||
OB[i][3] = opticalbench[i][2]
|
OB[i][3] = opticalbench[i][2]
|
||||||
|
|
||||||
elif opticalbench[i][3] == 'r':
|
elif opticalbench[i][3] == 'r':
|
||||||
if i == 0:
|
if i == 0:
|
||||||
position = position + opticalbench[i][2]
|
position = position + opticalbench[i][2]
|
||||||
OB[i][3] = position
|
OB[i][3] = position
|
||||||
|
|
||||||
else:
|
else:
|
||||||
if opticalbench[i-1][3] == 'a':
|
if opticalbench[i-1][3] == 'a':
|
||||||
position = position + opticalbench[i][2]
|
position = position + opticalbench[i][2]
|
||||||
OB[i][3] = position
|
OB[i][3] = position
|
||||||
elif opticalbench[i-1][3] == 'r':
|
elif opticalbench[i-1][3] == 'r':
|
||||||
position = position + opticalbench[i][2]
|
position = position + opticalbench[i][2]
|
||||||
OB[i][3] = position
|
OB[i][3] = position
|
||||||
elif opticalbench[i-1][3] == 'f':
|
elif opticalbench[i-1][3] == 'f':
|
||||||
OB[i][3] = position
|
OB[i][3] = position
|
||||||
else:
|
else:
|
||||||
print('Warning: Unknown position parameter.')
|
print('Warning: Unknown position parameter.')
|
||||||
print('Optical bench entry: '+str(i-1))
|
print('Optical bench entry: '+str(i-1))
|
||||||
break
|
break
|
||||||
else:
|
else:
|
||||||
print('Warning: Unknown position parameter.')
|
print('Warning: Unknown position parameter.')
|
||||||
print('Optical bench entry: '+str(i))
|
print('Optical bench entry: '+str(i))
|
||||||
break
|
break
|
||||||
|
|
||||||
|
|
||||||
#create ray-matrix
|
#create ray-matrix
|
||||||
M = [[0,0,0] for j in range(elements)]
|
M = [[0,0,0] for j in range(elements)]
|
||||||
M_tmp = np.empty((2,2))
|
M_tmp = np.empty((2,2))
|
||||||
M_tmp[:] = 0
|
M_tmp[:] = 0
|
||||||
|
|
||||||
#[start position, stop position, matrix] parameter is always z
|
#[start position, stop position, matrix] parameter is always z
|
||||||
for i in range(0,elements):
|
for i in range(0,elements):
|
||||||
M[i][0] = OB[i][2]
|
M[i][0] = OB[i][2]
|
||||||
M[i][1] = OB[i][3]
|
M[i][1] = OB[i][3]
|
||||||
|
|
||||||
if OB[i][0] == 'path':
|
if OB[i][0] == 'path':
|
||||||
if i == 0:
|
if i == 0:
|
||||||
M[i][2] = 'np.array([[1,(z-'+str(OB[i][2])+')/'+str(OB[i][1])+'],[0,1]])'
|
M[i][2] = 'np.array([[1,(z-'+str(OB[i][2])+')/'+str(OB[i][1])+'],[0,1]])'
|
||||||
else:
|
else:
|
||||||
z = OB[i][2]
|
z = OB[i][2]
|
||||||
M_tmp = eval(M[i-1][2])
|
M_tmp = eval(M[i-1][2])
|
||||||
M[i][2] = 'np.matmul(np.array([[1,(z-'+str(OB[i][2])+')/'+str(OB[i][1])+'],[0,1]]),\
|
M[i][2] = 'np.matmul(np.array([[1,(z-'+str(OB[i][2])+')/'+str(OB[i][1])+'],[0,1]]),\
|
||||||
np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']]))'
|
np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']]))'
|
||||||
|
|
||||||
|
|
||||||
elif OB[i][0] == 'thinlens':
|
elif OB[i][0] == 'thinlens':
|
||||||
if i == 0:
|
if i == 0:
|
||||||
M[i][2] = 'np.array([[1,0],[-1/'+str(OB[i][1])+',1]])'
|
M[i][2] = 'np.array([[1,0],[-1/'+str(OB[i][1])+',1]])'
|
||||||
|
|
||||||
else:
|
else:
|
||||||
z = OB[i][2]
|
z = OB[i][2]
|
||||||
M_tmp = eval('np.matmul(np.array([[1,0],[-1/'+str(OB[i][1])+',1]]),'+M[i-1][2]+')')
|
M_tmp = eval('np.matmul(np.array([[1,0],[-1/'+str(OB[i][1])+',1]]),'+M[i-1][2]+')')
|
||||||
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
||||||
|
|
||||||
|
|
||||||
elif OB[i][0] == 'thicklens':
|
elif OB[i][0] == 'thicklens':
|
||||||
# [1,0;n2-n1/(R2*n1),n2/n1] * [1,d/n20,1] * [1,0;(n1-n2)/(R1*n2),n1/n2]
|
# [1,0;n2-n1/(R2*n1),n2/n1] * [1,d/n20,1] * [1,0;(n1-n2)/(R1*n2),n1/n2]
|
||||||
M_tmp = np.matmul(np.array([[1,0],[(OB[i][1][1]-OB[i][1][0])/(OB[i][1][3]*OB[i][1][0]),OB[i][1][1]/OB[i][1][0]]]),
|
M_tmp = np.matmul(np.array([[1,0],[(OB[i][1][1]-OB[i][1][0])/(OB[i][1][3]*OB[i][1][0]),OB[i][1][1]/OB[i][1][0]]]),
|
||||||
np.matmul(np.array([[1,OB[i][1][4]/OB[i][1][1]],[0,1]]),
|
np.matmul(np.array([[1,OB[i][1][4]/OB[i][1][1]],[0,1]]),
|
||||||
np.array([[1,0],[(OB[i][1][0]-OB[i][1][1])/(OB[i][1][2]*OB[i][1][1]),OB[i][1][0]/OB[i][1][1]]])))
|
np.array([[1,0],[(OB[i][1][0]-OB[i][1][1])/(OB[i][1][2]*OB[i][1][1]),OB[i][1][0]/OB[i][1][1]]])))
|
||||||
if i == 1:
|
if i == 1:
|
||||||
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
||||||
else:
|
else:
|
||||||
z = OB[i][2]
|
z = OB[i][2]
|
||||||
M_tmp = eval('np.matmul(np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']]),'+M[i-1][2]+')')
|
M_tmp = eval('np.matmul(np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']]),'+M[i-1][2]+')')
|
||||||
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
||||||
|
|
||||||
elif OB[i][0] == 'flatmirror':
|
elif OB[i][0] == 'flatmirror':
|
||||||
M_tmp = np.array([[1,0],[0,1]])
|
M_tmp = np.array([[1,0],[0,1]])
|
||||||
if i == 1:
|
if i == 1:
|
||||||
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
||||||
else:
|
else:
|
||||||
z = OB[i][2]
|
z = OB[i][2]
|
||||||
M_tmp = eval('np.matmul(np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']]),'+M[i-1][2]+')')
|
M_tmp = eval('np.matmul(np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']]),'+M[i-1][2]+')')
|
||||||
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
||||||
|
|
||||||
elif OB[i][0] == 'curvedmirror':
|
elif OB[i][0] == 'curvedmirror':
|
||||||
if OB[i][1][2] == 0: #0=tangential, 1=sagittal
|
if OB[i][1][2] == 0: #0=tangential, 1=sagittal
|
||||||
M_tmp = np.array([[1,0],[-2/(OB[i][1][0]*np.cos(np.pi*OB[i][1][1]/180)),1]])
|
M_tmp = np.array([[1,0],[-2/(OB[i][1][0]*np.cos(np.pi*OB[i][1][1]/180)),1]])
|
||||||
else:
|
else:
|
||||||
M_tmp = np.array([[1,0],[-2*np.cos(np.pi*OB[i][1][1]/180)/OB[i][1][0],1]])
|
M_tmp = np.array([[1,0],[-2*np.cos(np.pi*OB[i][1][1]/180)/OB[i][1][0],1]])
|
||||||
|
|
||||||
|
|
||||||
if i == 1:
|
if i == 1:
|
||||||
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
||||||
else:
|
else:
|
||||||
z = OB[i][2]
|
z = OB[i][2]
|
||||||
M_tmp = eval('np.matmul(np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']]),'+M[i-1][2]+')')
|
M_tmp = eval('np.matmul(np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']]),'+M[i-1][2]+')')
|
||||||
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
||||||
|
|
||||||
elif OB[i][0] == 'flatrefraction':
|
elif OB[i][0] == 'flatrefraction':
|
||||||
#M_tmp = [1,00,OB[i,2](2)/OB[i,2](1)]
|
#M_tmp = [1,00,OB[i,2](2)/OB[i,2](1)]
|
||||||
M_tmp = np.array([[1,0],[0,OB[i][1][0]/OB[i][1][1]]])
|
M_tmp = np.array([[1,0],[0,OB[i][1][0]/OB[i][1][1]]])
|
||||||
if i == 1:
|
if i == 1:
|
||||||
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
||||||
else:
|
else:
|
||||||
z = OB[i][2]
|
z = OB[i][2]
|
||||||
M_tmp = eval('np.matmul(np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']]),'+M[i-1][2]+')')
|
M_tmp = eval('np.matmul(np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']]),'+M[i-1][2]+')')
|
||||||
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
||||||
|
|
||||||
elif OB[i][0] == 'curvedrefraction':
|
elif OB[i][0] == 'curvedrefraction':
|
||||||
# [1 , 0 ]
|
# [1 , 0 ]
|
||||||
# [(n1-n2)/R*n2 , n1/n2]
|
# [(n1-n2)/R*n2 , n1/n2]
|
||||||
M_tmp = np.array([[1,0],[(OB[i][1][0]-OB[i][1][1])/(OB[i][1][2]*OB[i][1][1]),OB[i][1][0]/OB[i][1][1]]])
|
M_tmp = np.array([[1,0],[(OB[i][1][0]-OB[i][1][1])/(OB[i][1][2]*OB[i][1][1]),OB[i][1][0]/OB[i][1][1]]])
|
||||||
if i == 1:
|
if i == 1:
|
||||||
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
||||||
else:
|
else:
|
||||||
z = OB[i][2]
|
z = OB[i][2]
|
||||||
M_tmp = eval('np.matmul(np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']]),'+M[i-1][2]+')')
|
M_tmp = eval('np.matmul(np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']]),'+M[i-1][2]+')')
|
||||||
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
M[i][2] = 'np.array([['+str(M_tmp[0][0])+','+str(M_tmp[0][1])+'],['+str(M_tmp[1][0])+','+str(M_tmp[1][1])+']])'
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
x = np.linspace(M[0][0],M[-1][1],N)
|
x = np.linspace(M[0][0],M[-1][1],N)
|
||||||
R = np.zeros((N,2))
|
R = np.zeros((N,2))
|
||||||
w = np.zeros((N,1))
|
w = np.zeros((N,1))
|
||||||
|
|
||||||
j = 0
|
j = 0
|
||||||
for i in range(np.size(x)):
|
for i in range(np.size(x)):
|
||||||
z = x[i]
|
z = x[i]
|
||||||
while M[j][1] < z:
|
while M[j][1] < z:
|
||||||
j = j+1
|
j = j+1
|
||||||
|
|
||||||
M_tmp = eval(M[j][2])
|
M_tmp = eval(M[j][2])
|
||||||
|
|
||||||
q = ( M_tmp[0][0] * q0 + M_tmp[0][1] ) / ( M_tmp[1][0] * q0 + M_tmp[1][1])
|
q = ( M_tmp[0][0] * q0 + M_tmp[0][1] ) / ( M_tmp[1][0] * q0 + M_tmp[1][1])
|
||||||
if np.real(q) == 0:
|
if np.real(q) == 0:
|
||||||
R[i] = 0
|
R[i] = 0
|
||||||
else:
|
else:
|
||||||
R[i] = np.real( q ) * (1 + ( np.imag( q ) / np.real( q ))**2)
|
R[i] = np.real( q ) * (1 + ( np.imag( q ) / np.real( q ))**2)
|
||||||
|
|
||||||
if OB[j][0] == 'path':
|
if OB[j][0] == 'path':
|
||||||
w[i] = np.sqrt(lam*np.abs(np.imag(q))/np.pi)*np.sqrt(1+(np.real(q)/np.imag(q))**2)/np.sqrt(OB[j][1])
|
w[i] = np.sqrt(lam*np.abs(np.imag(q))/np.pi)*np.sqrt(1+(np.real(q)/np.imag(q))**2)/np.sqrt(OB[j][1])
|
||||||
else:
|
else:
|
||||||
w[i] = np.sqrt(lam*np.abs(np.imag(q))/np.pi)*np.sqrt(1+(np.real(q)/np.imag(q))**2)
|
w[i] = np.sqrt(lam*np.abs(np.imag(q))/np.pi)*np.sqrt(1+(np.real(q)/np.imag(q))**2)
|
||||||
|
|
||||||
|
|
||||||
M_tot = M_tmp
|
M_tot = M_tmp
|
||||||
|
|
||||||
return [x,w,R,OB,M_tot]
|
return [x,w,R,OB,M_tot]
|
||||||
|
|
||||||
|
|
||||||
|
|||||||
Loading…
x
Reference in New Issue
Block a user