COMPASS  5.4.4
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  Args:
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  Args:
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  The arguments <pitch> and <x>, <y> must be in the same unit.
146  This unit shall be [pixels] if one expects that the returned value
147  <smallsize> is also in [pixels].
148 
149  Args:
150 
151  pitch: (float) : pitch of the DM expressed in pixels
152 
153  coupling: (float) : coupling of the actuators
154 
155  x: indices of influence function in relative position x local coordinates (float). 0 = top of the influence function
156 
157  y: indices of influence function in relative position y local coordinates (float). 0 = top of the influence function
158 
159  :return:
160 
161  influ: (np.ndarray(dims=3,dtype=np.float64)) : cube of the IF for each actuator
162  """
163  irc = 1.16136 + 2.97422 * coupling + \
164  (-13.2381) * coupling**2 + 20.4395 * coupling**3
165 
166  p1 = 4.49469 + (7.25509 + (-32.1948 + 17.9493 * coupling) * coupling) * coupling
167  p2 = 2.49456 + (-0.65952 + (8.78886 - 6.23701 * coupling) * coupling) * coupling
168 
169  tmp_c = 1.0 / np.abs(irc)
170  ccc = (coupling - 1. + tmp_c**p1) / (np.log(tmp_c) * tmp_c**p2)
171 
172  smallsize = int(np.ceil(2 * irc * pitch + 10))
173  if (smallsize % 2 != 0):
174  smallsize += 1
175  # clip
176  if (x is None or y is None):
177  return smallsize
178  else:
179  # normalized coordiantes in local ref frame
180  x = np.abs(x) / (irc * pitch)
181  y = np.abs(y) / (irc * pitch)
182 
183  x[x < 1e-8] = 1e-8
184  x[x > 2] = 2.
185  y[y < 1e-8] = 1e-8
186  y[y > 2] = 2.
187  tmp = (1. - x**p1 + ccc * np.log(x) * x**p2) * \
188  (1. - y**p1 + ccc * np.log(y) * y**p2)
189  tmp = tmp * (x <= 1.0) * (y <= 1.0)
190  return tmp
191 
192 
193 def makeRadialSchwartz(pitch: float, coupling: float, x=None, y=None):
194  """ Compute radial Schwartz influence function
195 
196  Args:
197 
198  pitch: (float) : pitch of the DM expressed in pixels
199 
200  coupling: (float) : coupling of the actuators
201 
202  x: indices of influence function in relative position x local coordinates (float). 0 = top of the influence function
203 
204  y: indices of influence function in relative position y local coordinates (float). 0 = top of the influence function
205 
206  :return:
207 
208  influ: (np.ndarray(dims=3,dtype=np.float64)) : cube of the IF for each actuator
209  """
210  k = 6 # order of the Schwartz function
211  #
212  a = pitch / np.sqrt(k / (np.log(coupling) - k) + 1.)
213  smallsize = int(2 * np.ceil(a) + 2)
214  if (x is None or y is None):
215  return smallsize
216  else:
217  r = (x * x + y * y) / (a * a)
218  ok = np.where(r < 1)
219  sc = np.zeros(r.shape)
220  sc[ok] = np.exp((k / (r[ok] - 1.0)) + k)
221  #influ[:,:,:] = sc[:,:,None] * np.ones(ntotact)[None,None,:]
222  return sc
223 
224 
225 def makeSquareSchwartz(pitch: float, coupling: float, x=None, y=None):
226  """ Compute Square Schwartz influence function
227 
228  Args:
229 
230  pitch: (float) : pitch of the DM expressed in pixels
231 
232  coupling: (float) : coupling of the actuators
233 
234  x: indices of influence function in relative position x local coordinates (float). 0 = top of the influence function
235 
236  y: indices of influence function in relative position y local coordinates (float). 0 = top of the influence function
237 
238  :return:
239 
240  influ: (np.ndarray(dims=3,dtype=np.float64)) : cube of the IF for each actuator
241  """
242  k = 6 # order of the Schwartz function
243  #
244  a = pitch / np.sqrt(k / (np.log(coupling) - k) + 1.)
245 
246  if (x is None or y is None):
247  smallsize = int(2 * np.ceil(a) + 2)
248  return smallsize
249  else:
250  xx = (x / a)**2
251  yy = (y / a)**2
252  ok = np.where((xx < 1) * (yy < 1))
253  sc = np.zeros(xx.shape)
254  sc[ok] = np.exp((k / (xx[ok] - 1)) + k) * \
255  np.exp((k / (yy[ok] - 1)) + k)
256  return sc
257 
258 
259 def makeBlacknutt(pitch: float, coupling: float, x=None, y=None):
260  """ Compute Blacknutt influence function
261  Attention, ici on ne peut pas choisir la valeur de coupling.
262  La variable a ete laissee dans le code juste pour compatibilité avec les
263  autres fonctions, mais elle n'est pas utilisee.
264 
265  Args:
266 
267  pitch: (float): pitch of the DM expressed in pixels
268 
269  coupling: (float) : coupling of the actuators
270 
271  x: indices of influence function in relative position x local coordinates (float). 0 = top of the influence function
272 
273  y: indices of influence function in relative position y local coordinates (float). 0 = top of the influence function
274 
275  :return:
276 
277  influ: (np.ndarray(dims=3,dtype=np.float64)) : cube of the IF for each actuator
278  """
279  smallsize = int(np.ceil(4 * pitch + 1))
280  if (x is None or y is None):
281  return smallsize
282  else:
283  cg = smallsize // 2
284  xx = x / float(cg)
285  yy = y / float(cg)
286  a = np.array([0.355768, 0.487396, 0.144232, 0.012604], dtype=np.float32)
287  ok = np.where((np.abs(xx) < 1) * (np.abs(yy) < 1))
288  sc = np.zeros(xx.shape)
289  sc[ok] = (a[0] + a[1] * np.cos(np.pi * xx[ok]) +
290  a[2] * np.cos(2 * np.pi * xx[ok]) + a[3] * np.cos(3 * np.pi * xx[ok])) *\
291  (a[0] + a[1] * np.cos(np.pi * yy[ok]) +
292  a[2] * np.cos(2 * np.pi * yy[ok]) + a[3] * np.cos(3 * np.pi * yy[ok]))
293 
294  return sc
295 
296 
297 def makeGaussian(pitch: float, coupling: float, x=None, y=None):
298  """ Compute Gaussian influence function. Coupling parameter is not taken into account
299 
300  Args:
301 
302  pitch: (float) : pitch of the DM expressed in pixels
303 
304  coupling: (float) : coupling of the actuators
305 
306  x: indices of influence function in relative position x local coordinates (float). 0 = top of the influence function
307 
308  y: indices of influence function in relative position y local coordinates (float). 0 = top of the influence function
309 
310  :return:
311 
312  influ: (np.ndarray(dims=3,dtype=np.float64)) : cube of the IF for each actuator
313  """
314  irc = 1.16136 + 2.97422 * coupling + \
315  (-13.2381) * coupling**2 + 20.4395 * coupling**3
316  tmp_c = 1.0 / np.abs(irc)
317 
318  smallsize = int(np.ceil(2 * irc * pitch + 10))
319  if (smallsize % 2 != 0):
320  smallsize += 1
321 
322  if (x is None or y is None):
323  return smallsize
324  else:
325  sig = pitch / np.sqrt(-2 * np.log(coupling))
326  gauss = np.exp(-0.5 * ((x / sig)**2 + (y / sig)**2))
327  # gauss /= gauss.max()
328  # xdg = np.linspace(-1, 1, smallsize, dtype=np.float32)
329  # x = np.tile(xdg, (smallsize, 1))
330  # y = x.T
331  # sig = 0.8
332  # gauss = 1 / np.cos(np.exp(-(x**2 / sig + y**2 / sig))**2)
333  # # Force value at zero on array limits
334  # gauss -= gauss[gauss.shape[0] // 2].min()
335  # gauss[gauss < 0.] = 0
336  # gauss /= gauss.max() # Normalize
337  return gauss
338 
339 
340 def makeBessel(pitch: float, coupling: float, x: np.ndarray = None, y: np.ndarray = None,
341  patternType: bytes = PatternType.SQUARE):
342  """ Compute Bessel influence function
343 
344  Args:
345 
346  pitch: (float) : pitch of the DM expressed in pixels
347 
348  coupling: (float) : coupling of the actuators
349 
350  x: indices of influence function in relative position x local coordinates (float). 0 = top of the influence function
351 
352  y: indices of influence function in relative position y local coordinates (float). 0 = top of the influence function
353 
354  :return:
355 
356  influ: (np.ndarray(dims=3,dtype=np.float64)) : cube of the IF for each actuator
357  """
358  smallsize = int(np.ceil(pitch * 3.2))
359 
360  if (x is None or y is None):
361  return smallsize
362  else:
363  # size_pitch = smallsize/np.float32(p_dm._pitch) # size of influence fonction in pitch
364  # xdg= np.linspace(-1*(size_pitch/3.2),size_pitch/3.2, smallsize,dtype=np.float32)
365  # x = np.tile(xdg, (smallsize,1))
366  # y = x.T
367  influ_u = bessel_influence(x / (1.6 * pitch), y / (1.6 * pitch), patternType)
368  influ_u = influ_u * (influ_u >= 0)
369  return influ_u
Numerical constants for shesha and config enumerations for safe-typing.
Definition: constants.py:1
def makeRadialSchwartz(float pitch, float coupling, x=None, y=None)
Compute radial Schwartz influence function.
Definition: influ_util.py:209
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:278
def makeRigaut(float pitch, float coupling, x=None, y=None)
Compute 'Rigaut-like' influence function.
Definition: influ_util.py:162
def besel_orth(m, n, phi, r)
TODO docstring.
Definition: influ_util.py:60
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:357
def bessel_influence(xx, yy, type_i=PatternType.SQUARE)
TODO docstring.
Definition: influ_util.py:86
def makeSquareSchwartz(float pitch, float coupling, x=None, y=None)
Compute Square Schwartz influence function.
Definition: influ_util.py:241
def makeGaussian(float pitch, float coupling, x=None, y=None)
Compute Gaussian influence function.
Definition: influ_util.py:313