COMPASS  5.0.0
End-to-end AO simulation tool using GPU acceleration
influ_util.py
1 
37 
38 import numpy as np
39 import scipy.special as sp
40 
41 from shesha.constants import DmType, PatternType
42 
43 
44 def besel_orth(m, n, phi, r):
45  """ TODO: docstring
46 
47  :parameters:
48 
49  m:
50 
51  n:
52 
53  phi:
54 
55  r:
56 
57  :return:
58 
59  B:
60  """
61  # fonction de bessel fourier orthogonale (BFOFS)
62  if (m == 0):
63  B = sp.jn(0, sp.jn_zeros(0, n)[n - 1] * r)
64  elif (m > 0):
65  B = sp.jn(m, sp.jn_zeros(m, n)[n - 1] * r) * np.sin(m * phi)
66  else:
67  B = sp.jn(np.abs(m),
68  sp.jn_zeros(np.abs(m), n)[n - 1] * r) * np.cos(np.abs(m) * phi)
69  return B
70 
71 
72 def bessel_influence(xx, yy, type_i=PatternType.SQUARE):
73  """ TODO: docstring
74 
75  :parameters:
76 
77  xx:
78 
79  yy:
80 
81  type_i: (optional)
82 
83  :return:
84 
85  influ
86  """
87 
88  # base sur l article numerical model of the influence function of deformable mirrors based on bessel Fourier orthogonal functions
89  # corespond a 3.2pitch
90 
91  influ = np.zeros(xx.shape, dtype=np.float32)
92 
93  # construction des tableaux :
94 
95  # construction des coordonnée cartesienne
96  # x = np.arange(size)-middle # -->
97  # y = (np.arange(size)-middle)*-1 # -->
98  #xx,yy = np.meshgrid(x,y)
99  # passage en coordonnée polaire
100  r = np.sqrt(xx**2 + yy**2)
101  phi = np.arctan2(yy, xx) # +(np.pi/8.) #petite correction de rotation
102 
103  # coef for square IF
104  a0 = [
105  0.3826, 0.5207, 0.2841, -0.0146, -0.1103, -0.0818, -0.0141, 0.0123, 0.0196,
106  0.0037
107  ]
108  am = [-0.0454, -0.1114, -0.1125, -0.0397, 0.0146, 0.0217, 0.0085, -0.0012, -0.0040]
109  a = [-0.0002, -0.0004, -0.0001, 0.0004, 0.0005, 0.0003, 0.0001, 0, 0]
110 
111  # search coef for hexagonal IF (m =0, 6, -6 --> 28 term)
112  # a0 ->10
113  # a6 ->9
114  # am6 ->9
115 
116  if type_i == PatternType.HEXA or type_i == PatternType.HEXAM4:
117  sym = 6
118 
119  else:
120  sym = 4
121 
122  # calcul pour m = 0
123  for i in range(len(a0)):
124  btemp = (a0[i] * besel_orth(0, i + 1, phi, r))
125 
126  influ = influ + btemp
127  #print("fin cas m=",0)
128 
129  # calcul pour m=6
130  for i in range(len(a)):
131  influ = influ + (a[i] * besel_orth(sym, i + 1, phi, r))
132  #print("fin cas m=",sym)
133 
134  # calcul pour m=-6
135  for i in range(len(am)):
136  influ = influ + (am[i] * besel_orth(-sym, i + 1, phi, r))
137  #print("fin cas m=",-sym)
138 
139  return influ
140 
141 
142 def makeRigaut(pitch: float, coupling: float, x=None, y=None):
143  """ Compute 'Rigaut-like' influence function
144 
145  :parameters:
146 
147  pitch: (float) : pitch of the DM expressed in pixels
148 
149  coupling: (float) : coupling of the actuators
150 
151  x: indices of influence function in relative position x local coordinates (float). 0 = top of the influence function
152 
153  y: indices of influence function in relative position y local coordinates (float). 0 = top of the influence function
154 
155  :return:
156 
157  influ: (np.ndarray(dims=3,dtype=np.float64)) : cube of the IF for each actuator
158  """
159  irc = 1.16136 + 2.97422 * coupling + \
160  (-13.2381) * coupling**2 + 20.4395 * coupling**3
161 
162  p1 = 4.49469 + (7.25509 + (-32.1948 + 17.9493 * coupling) * coupling) * coupling
163  p2 = 2.49456 + (-0.65952 + (8.78886 - 6.23701 * coupling) * coupling) * coupling
164 
165  tmp_c = 1.0 / np.abs(irc)
166  ccc = (coupling - 1. + tmp_c**p1) / (np.log(tmp_c) * tmp_c**p2)
167 
168  smallsize = int(np.ceil(2 * irc * pitch + 10))
169  if (smallsize % 2 != 0):
170  smallsize += 1
171  # clip
172  if (x is None or y is None):
173  return smallsize
174  else:
175  # normalized coordiantes in local ref frame
176  x = np.abs(x) / (irc * pitch)
177  y = np.abs(y) / (irc * pitch)
178 
179  x[x < 1e-8] = 1e-8
180  x[x > 2] = 2.
181  y[y < 1e-8] = 1e-8
182  y[y > 2] = 2.
183  tmp = (1. - x**p1 + ccc * np.log(x) * x**p2) * \
184  (1. - y**p1 + ccc * np.log(y) * y**p2)
185  tmp = tmp * (x <= 1.0) * (y <= 1.0)
186  return tmp
187 
188 
189 def makeRadialSchwartz(pitch: float, coupling: float, x=None, y=None):
190  """ Compute radial Schwartz influence function
191 
192  :parameters:
193 
194  pitch: (float) : pitch of the DM expressed in pixels
195 
196  coupling: (float) : coupling of the actuators
197 
198  x: indices of influence function in relative position x local coordinates (float). 0 = top of the influence function
199 
200  y: indices of influence function in relative position y local coordinates (float). 0 = top of the influence function
201 
202  :return:
203 
204  influ: (np.ndarray(dims=3,dtype=np.float64)) : cube of the IF for each actuator
205  """
206  k = 6 # order of the Schwartz function
207  #
208  a = pitch / np.sqrt(k / (np.log(coupling) - k) + 1.)
209  smallsize = int(2 * np.ceil(a) + 2)
210  if (x is None or y is None):
211  return smallsize
212  else:
213  r = (x * x + y * y) / (a * a)
214  ok = np.where(r < 1)
215  sc = np.zeros(r.shape)
216  sc[ok] = np.exp((k / (r[ok] - 1.0)) + k)
217  #influ[:,:,:] = sc[:,:,None] * np.ones(ntotact)[None,None,:]
218  return sc
219 
220 
221 def makeSquareSchwartz(pitch: float, coupling: float, x=None, y=None):
222  """ Compute Square Schwartz influence function
223 
224  :parameters:
225 
226  pitch: (float) : pitch of the DM expressed in pixels
227 
228  coupling: (float) : coupling of the actuators
229 
230  x: indices of influence function in relative position x local coordinates (float). 0 = top of the influence function
231 
232  y: indices of influence function in relative position y local coordinates (float). 0 = top of the influence function
233 
234  :return:
235 
236  influ: (np.ndarray(dims=3,dtype=np.float64)) : cube of the IF for each actuator
237  """
238  k = 6 # order of the Schwartz function
239  #
240  a = pitch / np.sqrt(k / (np.log(coupling) - k) + 1.)
241 
242  if (x is None or y is None):
243  smallsize = int(2 * np.ceil(a) + 2)
244  return smallsize
245  else:
246  xx = (x / a)**2
247  yy = (y / a)**2
248  ok = np.where((xx < 1) * (yy < 1))
249  sc = np.zeros(xx.shape)
250  sc[ok] = np.exp((k / (xx[ok] - 1)) + k) * \
251  np.exp((k / (yy[ok] - 1)) + k)
252  return sc
253 
254 
255 def makeBlacknutt(pitch: float, coupling: float, x=None, y=None):
256  """ Compute Blacknutt influence function
257  Attention, ici on ne peut pas choisir la valeur de coupling.
258  La variable a ete laissee dans le code juste pour compatibilité avec les
259  autres fonctions, mais elle n'est pas utilisee.
260 
261  :parameters:
262 
263  pitch: (float): pitch of the DM expressed in pixels
264 
265  coupling: (float) : coupling of the actuators
266 
267  x: indices of influence function in relative position x local coordinates (float). 0 = top of the influence function
268 
269  y: indices of influence function in relative position y local coordinates (float). 0 = top of the influence function
270 
271  :return:
272 
273  influ: (np.ndarray(dims=3,dtype=np.float64)) : cube of the IF for each actuator
274  """
275  smallsize = int(np.ceil(4 * pitch + 1))
276  if (x is None or y is None):
277  return smallsize
278  else:
279  cg = smallsize // 2
280  xx = x / float(cg)
281  yy = y / float(cg)
282  a = np.array([0.355768, 0.487396, 0.144232, 0.012604], dtype=np.float32)
283  ok = np.where((np.abs(xx) < 1) * (np.abs(yy) < 1))
284  sc = np.zeros(xx.shape)
285  sc[ok] = (a[0] + a[1] * np.cos(np.pi * xx[ok]) +
286  a[2] * np.cos(2 * np.pi * xx[ok]) + a[3] * np.cos(3 * np.pi * xx[ok])) *\
287  (a[0] + a[1] * np.cos(np.pi * yy[ok]) +
288  a[2] * np.cos(2 * np.pi * yy[ok]) + a[3] * np.cos(3 * np.pi * yy[ok]))
289 
290  return sc
291 
292 
293 def makeGaussian(pitch: float, coupling: float, x=None, y=None):
294  """ Compute Gaussian influence function. Coupling parameter is not taken into account
295 
296  :parameters:
297 
298  pitch: (float) : pitch of the DM expressed in pixels
299 
300  coupling: (float) : coupling of the actuators
301 
302  x: indices of influence function in relative position x local coordinates (float). 0 = top of the influence function
303 
304  y: indices of influence function in relative position y local coordinates (float). 0 = top of the influence function
305 
306  :return:
307 
308  influ: (np.ndarray(dims=3,dtype=np.float64)) : cube of the IF for each actuator
309  """
310  irc = 1.16136 + 2.97422 * coupling + \
311  (-13.2381) * coupling**2 + 20.4395 * coupling**3
312  tmp_c = 1.0 / np.abs(irc)
313 
314  smallsize = int(np.ceil(2 * irc * pitch + 10))
315  if (smallsize % 2 != 0):
316  smallsize += 1
317 
318  if (x is None or y is None):
319  return smallsize
320  else:
321  sig = pitch / np.sqrt(-2 * np.log(coupling))
322  gauss = np.exp(-0.5 * ((x / sig)**2 + (y / sig)**2))
323  # gauss /= gauss.max()
324  # xdg = np.linspace(-1, 1, smallsize, dtype=np.float32)
325  # x = np.tile(xdg, (smallsize, 1))
326  # y = x.T
327  # sig = 0.8
328  # gauss = 1 / np.cos(np.exp(-(x**2 / sig + y**2 / sig))**2)
329  # # Force value at zero on array limits
330  # gauss -= gauss[gauss.shape[0] // 2].min()
331  # gauss[gauss < 0.] = 0
332  # gauss /= gauss.max() # Normalize
333  return gauss
334 
335 
336 def makeBessel(pitch: float, coupling: float, x: np.ndarray = None, y: np.ndarray = None,
337  patternType: bytes = PatternType.SQUARE):
338  """ Compute Bessel influence function
339 
340  :parameters:
341 
342  pitch: (float) : pitch of the DM expressed in pixels
343 
344  coupling: (float) : coupling of the actuators
345 
346  x: indices of influence function in relative position x local coordinates (float). 0 = top of the influence function
347 
348  y: indices of influence function in relative position y local coordinates (float). 0 = top of the influence function
349 
350  :return:
351 
352  influ: (np.ndarray(dims=3,dtype=np.float64)) : cube of the IF for each actuator
353  """
354  smallsize = int(np.ceil(pitch * 3.2))
355 
356  if (x is None or y is None):
357  return smallsize
358  else:
359  # size_pitch = smallsize/np.float32(p_dm._pitch) # size of influence fonction in pitch
360  # xdg= np.linspace(-1*(size_pitch/3.2),size_pitch/3.2, smallsize,dtype=np.float32)
361  # x = np.tile(xdg, (smallsize,1))
362  # y = x.T
363  influ_u = bessel_influence(x / (1.6 * pitch), y / (1.6 * pitch), patternType)
364  influ_u = influ_u * (influ_u >= 0)
365  return influ_u
shesha.util.influ_util.besel_orth
def besel_orth(m, n, phi, r)
TODO docstring.
Definition: influ_util.py:60
shesha.constants
Numerical constants for shesha and config enumerations for safe-typing.
Definition: constants.py:1
shesha.util.influ_util.makeSquareSchwartz
def makeSquareSchwartz(float pitch, float coupling, x=None, y=None)
Compute Square Schwartz influence function.
Definition: influ_util.py:237
shesha.util.influ_util.makeBessel
def makeBessel(float pitch, float coupling, np.ndarray x=None, np.ndarray y=None, bytes patternType=PatternType.SQUARE)
Compute Bessel influence function.
Definition: influ_util.py:352
shesha.util.influ_util.makeGaussian
def makeGaussian(float pitch, float coupling, x=None, y=None)
Compute Gaussian influence function.
Definition: influ_util.py:309
shesha.util.influ_util.bessel_influence
def bessel_influence(xx, yy, type_i=PatternType.SQUARE)
TODO docstring.
Definition: influ_util.py:86
shesha.util.influ_util.makeRigaut
def makeRigaut(float pitch, float coupling, x=None, y=None)
Compute 'Rigaut-like' influence function.
Definition: influ_util.py:158
shesha.util.influ_util.makeBlacknutt
def makeBlacknutt(float pitch, float coupling, x=None, y=None)
Compute Blacknutt influence function Attention, ici on ne peut pas choisir la valeur de coupling.
Definition: influ_util.py:274
shesha.util.influ_util.makeRadialSchwartz
def makeRadialSchwartz(float pitch, float coupling, x=None, y=None)
Compute radial Schwartz influence function.
Definition: influ_util.py:205