COMPASS  5.4.4
End-to-end AO simulation tool using GPU acceleration
make_pupil.py
1 
37 
38 import numpy as np
39 import os
40 import scipy.ndimage.interpolation as interp
41 
42 from . import hdf5_util as h5u
43 from . import utilities as util
44 
45 from shesha.constants import ApertureType, SpiderType
46 
47 EELT_data = os.environ.get('SHESHA_ROOT') + "/data/apertures/"
48 
49 
50 def make_pupil(dim, pupd, tel, xc=-1, yc=-1, real=0, halfSpider=False):
51  """Initialize the system pupil
52 
53  Args:
54 
55  dim: (long) : = p_geom.pupdiam
56 
57  pupd: (long) : linear size of total pupil = p_geom.pupdiam
58 
59  tel: (Param_tel) : Telescope structure
60 
61  xc: (int) = p_geom.pupdiam / 2. - 0.5
62 
63  yc: (int) = p_geom.pupdiam / 2. - 0.5
64 
65  real: (int)
66 
67  TODO: complete
68  """
69 
70  if tel.type_ap == ApertureType.EELT_NOMINAL:
71  N_seg = 798
72  return make_EELT(dim, pupd, tel, N_seg)
73  elif (tel.type_ap == ApertureType.EELT):
74  return generateEeltPupilMask(dim, tel.t_spiders, xc, yc, tel.diam / dim, tel.gap,
75  tel.pupangle, D=tel.diam, halfSpider=halfSpider,
76  pitch=1.244683637214, nseg=33, inner_rad=4.1,
77  outer_rad=15.4, R=95.7853, nominalD=40,
78  half_seg=0.75, refl=tel.referr, nmissing=tel.nbrmissing)
79  elif (tel.type_ap == ApertureType.KECK):
80  seg_corner = 1.8
81  kpitch = seg_corner / 2 * np.sqrt(3)
82  knseg = 7
83  kinner_rad = 0.9
84  kouter_rad = 3.4
85  kR = 85
86  knominalD = 10.96
87  khalf_seg = 0.9
88  return generateEeltPupilMask(
89  dim, tel.t_spiders, xc, yc, tel.diam / dim, tel.gap, tel.pupangle,
90  D=tel.diam, cobs=tel.cobs, halfSpider=halfSpider, pitch=kpitch,
91  nseg=knseg, inner_rad=0.9, outer_rad=3.4, R=kR, nominalD=knominalD,
92  half_seg=0.9, refl=tel.referr, rotSpiderDegree=-30)
93  elif tel.type_ap == ApertureType.EELT_BP1:
94  print("ELT_pup_cobs = %5.3f" % 0.339)
95  N_seg = 768
96  return make_EELT(dim, pupd, tel, N_seg)
97  elif tel.type_ap == ApertureType.EELT_BP3:
98  print("ELT_pup_cobs = %5.3f" % 0.503)
99  N_seg = 672
100  return make_EELT(dim, pupd, tel, N_seg)
101  elif tel.type_ap == ApertureType.EELT_BP5:
102  print("ELT_pup_cobs = %5.3f" % 0.632)
103  N_seg = 558
104  return make_EELT(dim, pupd, tel, N_seg)
105  elif tel.type_ap == ApertureType.EELT_CUSTOM:
106  print("EELT_CUSTOM not implemented. Falling back to EELT-Nominal")
107  tel.set_type_ap(ApertureType.EELT_NOMINAL)
108  return make_EELT(dim, pupd, tel, N_seg)
109  elif tel.type_ap == ApertureType.VLT:
110  tel.set_cobs(0.14)
111  print("force_VLT_pup_cobs = %5.3f" % 0.14)
112  return make_VLT(dim, pupd, tel)
113  elif tel.type_ap == ApertureType.VLT_NOOBS:
114  tel.set_cobs(0.)
115  print("Remove central obstruction to the VLT")
116  return make_VLT(dim, pupd, tel)
117  elif tel.type_ap == ApertureType.GENERIC:
118  return make_pupil_generic(dim, pupd, tel.t_spiders, tel.spiders_type, xc, yc,
119  real, tel.cobs)
120  else:
121  raise NotImplementedError("Missing Pupil type.")
122 
123 
124 def make_pupil_generic(dim, pupd, t_spiders=0.01, spiders_type=SpiderType.SIX, xc=0,
125  yc=0, real=0, cobs=0):
126  """
127  Initialize the system pupil
128 
129  Args:
130 
131  dim: (long) : linear size of ???
132 
133  pupd: (long) : linear size of total pupil
134 
135  t_spiders: (float) : secondary supports ratio.
136 
137  spiders_type: (str) : secondary supports type: "four" or "six".
138 
139  xc: (int)
140 
141  yc: (int)
142 
143  real: (int)
144 
145  cobs: (float) : central obstruction ratio.
146 
147  TODO: complete
148  """
149 
150  pup = util.dist(dim, xc, yc)
151 
152  if (real == 1):
153  pup = np.exp(-(pup / (pupd * 0.5))**60.0)**0.69314
154  else:
155  pup = (pup < (pupd + 1.) / 2).astype(np.float32)
156 
157  if (cobs > 0):
158  if (real == 1):
159  pup -= np.exp(-(util.dist(dim, xc, yc) / (pupd * cobs * 0.5))**60.)**0.69314
160  else:
161  pup -= (util.dist(dim, xc, yc) < (pupd * cobs + 1.) * 0.5).astype(np.float32)
162 
163  step = 1. / dim
164  first = 0.5 * (step - 1)
165 
166  X = np.tile(np.arange(dim, dtype=np.float32) * step + first, (dim, 1))
167 
168  if (t_spiders < 0):
169  t_spiders = 0.01
170  t_spiders = t_spiders * pupd / dim
171 
172  if (spiders_type == "four"):
173 
174  s4_2 = 2 * np.sin(np.pi / 4)
175  t4 = np.tan(np.pi / 4)
176 
177  spiders_map = ((X.T > (X + t_spiders / s4_2) * t4) +
178  (X.T < (X - t_spiders / s4_2) * t4)).astype(np.float32)
179  spiders_map *= ((X.T > (-X + t_spiders / s4_2) * t4) +
180  (X.T < (-X - t_spiders / s4_2) * t4)).astype(np.float32)
181 
182  pup = pup * spiders_map
183 
184  elif (spiders_type == "six"):
185 
186  #angle = np.pi/(180/15.)
187  angle = 0
188  s2ma_2 = 2 * np.sin(np.pi / 2 - angle)
189  s6pa_2 = 2 * np.sin(np.pi / 6 + angle)
190  s6ma_2 = 2 * np.sin(np.pi / 6 - angle)
191  t2ma = np.tan(np.pi / 2 - angle)
192  t6pa = np.tan(np.pi / 6 + angle)
193  t6ma = np.tan(np.pi / 6 - angle)
194 
195  spiders_map = ((X.T > (-X + t_spiders / s2ma_2) * t2ma) +
196  (X.T < (-X - t_spiders / s2ma_2) * t2ma))
197  spiders_map *= ((X.T > (X + t_spiders / s6pa_2) * t6pa) +
198  (X.T < (X - t_spiders / s6pa_2) * t6pa))
199  spiders_map *= ((X.T > (-X + t_spiders / s6ma_2) * t6ma) +
200  (X.T < (-X - t_spiders / s6ma_2) * t6ma))
201  pup = pup * spiders_map
202 
203  print("Generic pupil created")
204  return pup
205 
206 
207 def make_VLT(dim, pupd, tel):
208  """
209  Initialize the VLT pupil
210 
211  Args:
212 
213  dim: (long) : linear size of ???
214 
215  pupd: (long) : linear size of total pupil
216 
217  tel: (Param_tel) : Telescope structure
218  """
219 
220  if (tel.t_spiders == -1):
221  print("force t_spider =%5.3f" % (0.09 / 18.))
222  tel.set_t_spiders(0.09 / 18.)
223  angle = 50.5 * np.pi / 180. # --> 50.5 degre *2 d'angle entre les spiders
224 
225  Range = (0.5 * (1) - 0.25 / dim)
226  X = np.tile(np.linspace(-Range, Range, dim, dtype=np.float32), (dim, 1))
227 
228  R = np.sqrt(X**2 + (X.T)**2)
229 
230  pup = ((R < 0.5) & (R > (tel.cobs / 2))).astype(np.float32)
231 
232  if (tel.t_spiders == -1):
233  print('No spider')
234  else:
235  spiders_map = (
236  (X.T >
237  (X - tel.cobs / 2 + tel.t_spiders / np.sin(angle)) * np.tan(angle)) +
238  (X.T < (X - tel.cobs / 2) * np.tan(angle))) * (X > 0) * (X.T > 0)
239  spiders_map += np.fliplr(spiders_map)
240  spiders_map += np.flipud(spiders_map)
241  spiders_map = interp.rotate(spiders_map, tel.pupangle, order=0, reshape=False)
242 
243  pup = pup * spiders_map
244 
245  print("VLT pupil created")
246  return pup
247 
248 
249 def make_EELT(dim, pupd, tel, N_seg=-1):
250  """
251  Initialize the EELT pupil
252 
253  Args:
254 
255  dim: (long) : linear size of ???
256 
257  pupd: (long) : linear size of total pupil
258 
259  tel: (Param_tel) : Telescope structure
260 
261  N_seg: (int)
262 
263  TODO: complete
264  TODO : add force rescal pup elt
265  """
266  if (N_seg == -1):
267  EELT_file = EELT_data + "EELT-Custom_N" + str(dim) + "_COBS" + str(
268  100 * tel.cobs) + "_CLOCKED" + str(tel.pupangle) + "_TSPIDERS" + str(
269  100 * tel.t_spiders) + "_MS" + str(
270  tel.nbrmissing) + "_REFERR" + str(
271  100 * tel.referr) + ".h5"
272  else:
273  EELT_file = EELT_data + tel.type_ap.decode('UTF-8') + "_N" + str(
274  dim) + "_COBS" + str(100 * tel.cobs) + "_CLOCKED" + str(
275  tel.pupangle) + "_TSPIDERS" + str(
276  100 * tel.t_spiders) + "_MS" + str(
277  tel.nbrmissing) + "_REFERR" + str(
278  100 * tel.referr) + ".h5"
279  if (os.path.isfile(EELT_file)):
280  print("reading EELT pupil from file ", EELT_file)
281  pup = h5u.readHdf5SingleDataset(EELT_file)
282  else:
283  print("creating EELT pupil...")
284  file = EELT_data + "Coord_" + tel.type_ap.decode('UTF-8') + ".dat"
285  data = np.fromfile(file, sep="\n")
286  data = np.reshape(data, (data.size // 2, 2))
287  x_seg = data[:, 0]
288  y_seg = data[:, 1]
289 
290  file = EELT_data + "EELT_MISSING_" + tel.type_ap.decode('UTF-8') + ".dat"
291  k_seg = np.fromfile(file, sep="\n").astype(np.int32)
292 
293  W = 1.45 * np.cos(np.pi / 6)
294 
295  #tel.set_diam(39.146)
296  #tel.set_diam(37.)
297  Range = (0.5 * (tel.diam * dim / pupd) - 0.25 / dim)
298  X = np.tile(np.linspace(-Range, Range, dim, dtype=np.float32), (dim, 1))
299 
300  if (tel.t_spiders == -1):
301  print("force t_spider =%5.3f" % (0.014))
302  tel.set_t_spiders(0.014)
303  #t_spiders=0.06
304  #tel.set_t_spiders(t_spiders)
305 
306  if (tel.nbrmissing > 0):
307  k_seg = np.sort(k_seg[:tel.nbrmissing])
308 
309  file = EELT_data + "EELT_REF_ERROR" + ".dat"
310  ref_err = np.fromfile(file, sep="\n")
311 
312  #mean_ref = np.sum(ref_err)/798.
313  #std_ref = np.sqrt(1./798.*np.sum((ref_err-mean_ref)**2))
314  #mean_ref=np.mean(ref_err)
315  std_ref = np.std(ref_err)
316 
317  ref_err = ref_err * tel.referr / std_ref
318 
319  if (tel.nbrmissing > 0):
320  ref_err[k_seg] = 1.0
321 
322  pup = np.zeros((dim, dim))
323 
324  t_3 = np.tan(np.pi / 3.)
325  if N_seg == -1:
326 
327  vect_seg = tel.vect_seg
328  for i in vect_seg:
329  Xt = X + x_seg[i]
330  Yt = X.T + y_seg[i]
331  pup+=(1-ref_err[i])*(Yt<0.5*W)*(Yt>=-0.5*W)*(0.5*(Yt+t_3*Xt)<0.5*W) \
332  *(0.5*(Yt+t_3*Xt)>=-0.5*W)*(0.5*(Yt-t_3*Xt)<0.5*W) \
333  *(0.5*(Yt-t_3*Xt)>=-0.5*W)
334 
335  else:
336  for i in range(N_seg):
337  Xt = X + x_seg[i]
338  Yt = X.T + y_seg[i]
339  pup+=(1-ref_err[i])*(Yt<0.5*W)*(Yt>=-0.5*W)*(0.5*(Yt+t_3*Xt)<0.5*W) \
340  *(0.5*(Yt+t_3*Xt)>=-0.5*W)*(0.5*(Yt-t_3*Xt)<0.5*W) \
341  *(0.5*(Yt-t_3*Xt)>=-0.5*W)
342  if (tel.t_spiders == 0):
343  print('No spider')
344  else:
345  t_spiders = tel.t_spiders * (tel.diam * dim / pupd)
346 
347  s2_6 = 2 * np.sin(np.pi / 6)
348  t_6 = np.tan(np.pi / 6)
349 
350  spiders_map = np.abs(X) > t_spiders / 2
351  spiders_map *= (
352  (X.T >
353  (X + t_spiders / s2_6) * t_6) + (X.T <
354  (X - t_spiders / s2_6) * t_6))
355  spiders_map *= (
356  (X.T >
357  (-X + t_spiders / s2_6) * t_6) + (X.T <
358  (-X - t_spiders / s2_6) * t_6))
359 
360  pup = pup * spiders_map
361 
362  if (tel.pupangle != 0):
363  pup = interp.rotate(pup, tel.pupangle, reshape=False, order=0)
364 
365  print("writing EELT pupil to file ", EELT_file)
366  h5u.writeHdf5SingleDataset(EELT_file, pup)
367 
368  print("EELT pupil created")
369  return pup
370 
371 
372 def make_phase_ab(dim, pupd, tel, pup=None, xc=-1, yc=-1, real=0, halfSpider=False):
373  """Compute the EELT M1 phase aberration
374 
375  Args:
376 
377  dim: (long) : linear size of ???
378 
379  pupd: (long) : linear size of total pupil
380 
381  tel: (Param_tel) : Telescope structure
382 
383  pup: (?)
384 
385  TODO: complete
386  """
387 
388  if ((tel.type_ap == ApertureType.GENERIC) or (tel.type_ap == ApertureType.VLT)):
389  return np.zeros((dim, dim)).astype(np.float32)
390 
391  elif tel.type_ap == ApertureType.EELT:
392 
393  return generateEeltPupilMask(
394  dim, 0, xc, yc, tel.diam / dim, tel.gap, tel.pupangle, D=tel.diam,
395  halfSpider=halfSpider, pitch=1.244683637214, nseg=33, inner_rad=4.1,
396  outer_rad=15.4, R=95.7853, nominalD=40, half_seg=0.75,
397  refl=[tel.std_piston, tel.std_tt, tel.std_tt])
398  elif (tel.type_ap == ApertureType.KECK):
399  seg_corner = 1.8
400  kpitch = seg_corner / 2 * np.sqrt(3)
401  knseg = 7
402  kinner_rad = 0.9
403  kouter_rad = 3.4
404  kR = 85
405  knominalD = 10.96
406  khalf_seg = 0.9
407  return generateEeltPupilMask(
408  dim, 0, xc, yc, tel.diam / dim, tel.gap, tel.pupangle, D=tel.diam,
409  cobs=tel.cobs, halfSpider=halfSpider, pitch=kpitch, nseg=knseg,
410  inner_rad=0.9, outer_rad=3.4, R=kR, nominalD=knominalD, half_seg=0.9,
411  refl=[tel.std_piston, tel.std_tt, tel.std_tt])
412  else:
413  ab_file = EELT_data + "aberration_" + tel.type_ap.decode('UTF-8') + \
414  "_N" + str(dim) + "_NPUP" + str(np.where(pup)[0].size) + "_CLOCKED" + str(
415  tel.pupangle) + "_TSPIDERS" + str(
416  100 * tel.t_spiders) + "_MS" + str(tel.nbrmissing) + "_REFERR" + str(
417  100 * tel.referr) + "_PIS" + str(
418  tel.std_piston) + "_TT" + str(tel.std_tt) + ".h5"
419  if (os.path.isfile(ab_file)):
420  print("reading aberration phase from file ", ab_file)
421  phase_error = h5u.readHdf5SingleDataset(ab_file)
422  else:
423  print("computing M1 phase aberration...")
424 
425  std_piston = tel.std_piston
426  std_tt = tel.std_tt
427 
428  W = 1.45 * np.cos(np.pi / 6)
429 
430  file = EELT_data + "EELT_Piston_" + tel.type_ap.decode('UTF-8') + ".dat"
431  p_seg = np.fromfile(file, sep="\n")
432  #mean_pis=np.mean(p_seg)
433  std_pis = np.std(p_seg)
434  p_seg = p_seg * std_piston / std_pis
435  N_seg = p_seg.size
436 
437  file = EELT_data + "EELT_TT_" + tel.type_ap.decode('UTF-8') + ".dat"
438  tt_seg = np.fromfile(file, sep="\n")
439 
440  file = EELT_data + "EELT_TT_DIRECTION_" + tel.type_ap.decode('UTF-8') + ".dat"
441  tt_phi_seg = np.fromfile(file, sep="\n")
442 
443  phase_error = np.zeros((dim, dim))
444  phase_tt = np.zeros((dim, dim))
445  phase_defoc = np.zeros((dim, dim))
446 
447  file = EELT_data + "Coord_" + tel.type_ap.decode('UTF-8') + ".dat"
448  data = np.fromfile(file, sep="\n")
449  data = np.reshape(data, (data.size // 2, 2))
450  x_seg = data[:, 0]
451  y_seg = data[:, 1]
452 
453  Range = (0.5 * (tel.diam * dim / pupd) - 0.25 / dim)
454 
455  X = np.tile(np.linspace(-Range, Range, dim, dtype=np.float32), (dim, 1))
456  t_3 = np.tan(np.pi / 3.)
457 
458  for i in range(N_seg):
459  Xt = X + x_seg[i]
460  Yt = X.T + y_seg[i]
461  SEG=(Yt<0.5*W)*(Yt>=-0.5*W)*(0.5*(Yt+t_3*Xt)<0.5*W) \
462  *(0.5*(Yt+t_3*Xt)>=-0.5*W)*(0.5*(Yt-t_3*Xt)<0.5*W) \
463  *(0.5*(Yt-t_3*Xt)>=-0.5*W)
464 
465  if (i == 0):
466  N_in_seg = np.sum(SEG)
467  Hex_diam = 2 * np.max(
468  np.sqrt(Xt[np.where(SEG)]**2 + Yt[np.where(SEG)]**2))
469 
470  if (tt_seg[i] != 0):
471  TT = tt_seg[i] * (
472  np.cos(tt_phi_seg[i]) * Xt + np.sin(tt_phi_seg[i]) * Yt)
473  mean_tt = np.sum(TT[np.where(SEG == 1)]) / N_in_seg
474  phase_tt += SEG * (TT - mean_tt)
475 
476  #TODO defocus
477 
478  phase_error += SEG * p_seg[i]
479 
480  N_EELT = np.where(pup)[0].size
481  if (np.sum(phase_tt) != 0):
482  phase_tt *= std_tt / np.sqrt(
483  1. / N_EELT * np.sum(phase_tt[np.where(pup)]**2))
484 
485  #TODO defocus
486 
487  phase_error += phase_tt + phase_defoc
488 
489  if (tel.pupangle != 0):
490  phase_error = interp.rotate(phase_error, tel.pupangle, reshape=False,
491  order=2)
492 
493  print("phase aberration created")
494  print("writing aberration filel to file ", ab_file)
495  h5u.writeHdf5SingleDataset(ab_file, phase_error)
496 
497  return phase_error
498 
499 
500 """
501 ooooooooooo ooooo ooooooooooo oooooooooo ooooo oooooooo8 ooooooo
502  888 88 888 88 888 88 888 888 888 o888 88 o888 888o
503  888ooo8 888 888 888oooo88 888 888 888 888
504  888 oo 888 o 888 888 88o 888 888o oo 888o o888
505 o888ooo8888 o888ooooo88 o888o o888o 88o8 o888o 888oooo88 88ooo88
506 """
507 
508 
509 def generateEeltPupilMask(npt, dspider, i0, j0, pixscale, gap, rotdegree, D=40.0, cobs=0,
510  centerMark=0, halfSpider=False, pitch=1.244683637214, nseg=33,
511  inner_rad=4.1, outer_rad=15.4, R=95.7853, nominalD=40,
512  half_seg=0.75, refl=None, rotSpiderDegree=None, nmissing=0):
513  """
514  Generates a boolean pupil mask of the binary EELT pupil
515  on a map of size (npt, npt).
516 
517 
518  :returns: pupil image (npt, npt), boolean
519  :param int npt: size of the output array
520  :param float dspider: width of spiders in meters
521  :param float i0, j0: index of pixels where the pupil should be centred = p_geom.pupdiam / 2. - 0.5
522  Can be floating-point indexes.
523  :param float pixscale: size of a pixel of the image, in meters = ptel.diam/(p_geom.pupdiam / 2. - 0.5)
524  :param float gap: gap between 2 segments in metres
525  :param float rotdegree: rotation angle of the pupil, in degrees.
526  :param float D: diameter of the pupil. For the nominal EELT, D shall
527  be set to 40.0
528  :param int centerMark: when centerMark!=0, a pixel is added at the centre of
529  symmetry of the pupil in order to debug things using compass.
530  centerMark==1 draws a point
531  centerMark==2 draws 2 lines
532  :param bool halfSpider: half Spider computation flag
533  :param float pitch: segment pitch
534  :param int nseg: number of segments across the diameter
535  :param float inner_rad: Inner radius [meters]
536  :param float outter_rad: outter radius [meters]
537  :param float R: M1 curvature radius
538  :param float nominalD: diameter needed to get nominal aperture after projection
539  :param float half_seg: segment half size
540  :param float refl: std of the reflectivity of each segment
541 
542  :Example:
543  npt = p_geom.pupdiam
544  D = p_tel.diam
545  i0 = npt / 2. - 0.5
546  j0 = npt / 2. - 0.5
547  rotdegree = 0.
548  pixscale = D/(npt / 2. - 0.5)
549  dspider = 0.51
550  gap = 0.0
551  pup = generateEeltPupilMask(npt, dspider, i0, j0, pixscale, gap, rotdegree)
552 
553  """
554  rot = rotdegree * np.pi / 180
555 
556  if rotSpiderDegree is None:
557  rotSpider = rot
558  else:
559  rotSpider = rotSpiderDegree * np.pi / 180
560 
561  # Generation of segments coordinates.
562  # hx and hy have a shape [6,798] describing the 6 vertex of the 798
563  # hexagonal mirrors
564  #hx, hy = generateCoordSegments( D, rot)
565  hx, hy = generateCoordSegments(D, rot, pitch=pitch, nseg=nseg, inner_rad=inner_rad,
566  outer_rad=outer_rad, R=R, nominalD=nominalD)
567  # From the data of hex mirrors, we build the pupil image using
568  # boolean
569  #pup = generateSegmentProperties(True, hx, hy, i0, j0, pixscale, gap, npt, D)
570  np.random.seed(42) #to ensure we have the same distribution for different simulations
571  if (refl == 0):
572  refl = True
573  elif np.isscalar(refl):
574  referr = np.random.random(hx.size)
575  referr = referr * refl / np.std(referr)
576  refl = np.ones(hx.size) - referr
577  elif type(refl) == list:
578  if len(refl) == 3:
579  refpist = np.random.randn(hx.size)
580  refpist = refpist * refl[0] / np.std(refpist)
581  reftip = np.random.randn(hx.size)
582  reftip = reftip * refl[1] / np.std(reftip)
583  reftilt = np.random.randn(hx.size)
584  reftilt = reftilt * refl[2] / np.std(reftilt)
585  refl = np.array([refpist, reftip, reftilt])
586  else:
587  raise ValueError(
588  "refl param must be None, scalar (reflectivity std error) or list of 3 elements (piston, tip and tilt std errors)"
589  )
590 
591  if (nmissing != 0):
592  np.random.seed(44)
593  ind_missing = np.random.randint(0, hx.size / 6, nmissing)
594  #ind_missing = np.arange(798)
595  print(ind_missing)
596  if np.isscalar(refl):
597  refl = np.ones(hx.size)
598  refl[ind_missing] = 0
599  #refl[ind_missing] = (ind_missing % 133 / (hx.size / 6)).astype(float)
600  elif (refl.shape[0] == 3): #piston tip tilt
601  print("WARNING: there are piston, tip and tilt std errors, missing segments will be ignored")
602  else: #there is already a 1D array containing reflectivity error
603  refl[ind_missing] = 0
604 
605 
606  pup = generateSegmentProperties(refl, hx, hy, i0, j0, pixscale, gap, npt, D,
607  nominalD=nominalD, pitch=pitch, half_seg=half_seg)
608  # SPIDERS ............................................
609  nspider = 3 # for the day where we have more/less spiders ;-)
610  if (dspider > 0 and nspider > 0):
611  if (halfSpider is True):
612  pup = pup * fillHalfSpider(npt, nspider, dspider, i0, j0, pixscale,
613  rotSpider)
614  else:
615  pup = pup * fillSpider(npt, nspider, dspider, i0, j0, pixscale, rotSpider)
616 
617  # Rajout d'un pixel au centre (pour marquer le centre) ou d'une croix,
618  # selon la valeur de centerMark
619  if centerMark:
620  pup = np.logical_xor(pup, centrePourVidal(npt, i0, j0, centerMark))
621 
622  if cobs > 0:
623  obstru = (util.dist(pup.shape[0], pup.shape[0] // 2 + 0.5,
624  pup.shape[0] // 2 + 0.5) >=
625  (pup.shape[0] * cobs + 1.) * 0.5).astype(np.float32)
626  pup = pup * obstru
627  return pup
628 
629 
630 def fillPolygon(x, y, i0, j0, scale, gap, N, index=0):
631  """
632  From a list of points defined by their 2 coordinates list
633  x and y, creates a filled polygon with sides joining the points.
634  The polygon is created in an image of size (N, N).
635  The origin (x,y)=(0,0) is mapped at pixel i0, j0 (both can be
636  floating-point values).
637  Arrays x and y are supposed to be in unit U, and scale is the
638  pixel size in U units.
639 
640  :returns: filled polygon (N, N), boolean
641  :param float x, y: list of points defining the polygon
642  :param float i0, j0: index of pixels where the pupil should be centred.
643  Can be floating-point indexes.
644  :param float scale: size of a pixel of the image, in same unit as x and y.
645  :param float N: size of output image.
646 
647  :Example:
648  x = np.array([1,-1,-1.5,0,1.1])
649  y = np.array([1,1.5,-0.2,-2,0])
650  N = 200
651  i0 = N/2
652  j0 = N/2
653  gap = 0.
654  scale = 0.03
655  pol = fillPolygon(x, y, i0, j0, scale, gap, N, index=2)
656 
657  """
658  # define coordinates map centred on (i0,j0) with same units as x,y.
659  X = (np.arange(N) - i0) * scale
660  Y = (np.arange(N) - j0) * scale
661  X, Y = np.meshgrid(X, Y, indexing='ij') # indexage [x,y]
662 
663  # define centre x0, y0 of the polygon
664  x0 = np.mean(x)
665  y0 = np.mean(y)
666 
667  # compute angles of all pixels coordinates of the map, and all
668  # corners of the polygon
669  T = (np.arctan2(Y - y0, X - x0) + 2 * np.pi) % (2 * np.pi)
670  t = (np.arctan2(y - y0, x - x0) + 2 * np.pi) % (2 * np.pi)
671 
672  # on va voir dans quel sens ca tourne. Je rajoute ca pour que ca marche
673  # quel que soit le sens de rotation des points du polygone.
674  # En fait, j'aurais peut etre pu classer les points par leur angle, pour
675  # etre sur que ca marche meme si les points sont donnes dans ts les cas
676  sens = np.median(np.diff(t))
677  if sens < 0:
678  x = x[::-1]
679  y = y[::-1]
680  t = t[::-1]
681 
682  # re-organise order of polygon points so that it starts from
683  # angle = 0, or at least closest to 0.
684  imin = t.argmin() # position of the minimum
685  if imin != 0:
686  x = np.roll(x, -imin)
687  y = np.roll(y, -imin)
688  t = np.roll(t, -imin)
689 
690  # For each couple of consecutive corners A, B, of the polygon, one fills
691  # the triangle AOB with True.
692  # Last triangle has a special treatment because it crosses the axis
693  # with theta=0=2pi
694  n = x.shape[0] # number of corners of polygon
695  indx, indy = (np.array([], dtype=np.int64), np.array([], dtype=np.int64))
696  distedge = np.array([], dtype=np.float32)
697  for i in range(n):
698  j = i + 1 # j=element next i except when i==n : then j=0 (cycling)
699  if j == n:
700  j = 0
701  sub = np.where((T >= t[-1]) | (T <= (t[0])))
702  else:
703  sub = np.where((T >= t[i]) & (T <= t[j]))
704  # compute unitary vector des 2 sommets
705  dy = y[j] - y[i]
706  dx = x[j] - x[i]
707  vnorm = np.sqrt(dx**2 + dy**2)
708  dx /= vnorm
709  dy /= vnorm
710  # calcul du produit vectoriel
711  crossprod = dx * (Y[sub] - y[i]) - dy * (X[sub] - x[i])
712  tmp = crossprod > gap
713  indx = np.append(indx, sub[0][tmp])
714  indy = np.append(indy, sub[1][tmp])
715  distedge = np.append(distedge, crossprod[tmp])
716 
717  # choice of what is returned : either only the indexes, or the
718  # boolean map
719  if index == 1:
720  return (indx, indy, distedge)
721  elif index == 2:
722  a = np.zeros((N, N))
723  a[indx, indy] = distedge
724  return a
725  else:
726  a = np.zeros((N, N), dtype=bool)
727  a[indx, indy] = True # convention [x,y]
728 
729  return a
730 
731 
732 def compute6Segments(pupNoSpiders, N, pixscale, dspider, i0, j0, rot=0):
733  """
734  N = p_geom.pupdiam
735  i0 = j0 = N / 2. - 0.5
736  D = p_tel.diam
737  pixscale = D/N
738  dspider = 0.51
739 
740  Utilisee dans compass/shesha/shesha/supervisor/canapassSupervisor.py
741  pour le slaving des actus.
742  """
743  npts = 200
744  msk = np.zeros((6, pupNoSpiders.shape[0], pupNoSpiders.shape[1]))
745  mid = int(pupNoSpiders.shape[0] / 2)
746  tmp = np.zeros((pupNoSpiders.shape[0], pupNoSpiders.shape[1]))
747  for angle in np.linspace(0, 60, npts):
748  tmp += compute1Spider(0, N, dspider, i0, j0, pixscale, angle * 2 * np.pi / 360)
749  msk[5, :, :] = tmp < npts * pupNoSpiders
750  msk[5, mid:, :] = 0
751  msk[2, :, :] = tmp < npts * pupNoSpiders
752  msk[2, :mid, :] = 0
753  tmp *= 0
754  for angle in np.linspace(0, 60, npts):
755  tmp += compute1Spider(1, N, dspider, i0, j0, pixscale, angle * 2 * np.pi / 360)
756  msk[0, :, :] = tmp < npts * pupNoSpiders
757  msk[0, mid:, :] = 0
758  msk[3, :, :] = tmp < npts * pupNoSpiders
759  msk[3, :mid, :] = 0
760  tmp *= 0
761  for angle in np.linspace(0, 60, npts):
762  tmp += compute1Spider(2, N, dspider, i0, j0, pixscale, angle * 2 * np.pi / 360)
763  msk[1, :, :] = tmp < npts * pupNoSpiders
764  msk[1, :, :mid] = 0
765  msk[4, :, :] = tmp < npts * pupNoSpiders
766  msk[4, :, mid:] = 0
767  return msk
768 
769 
770 def compute1Spider(nspider, N, dspider, i0, j0, scale, rot):
771  """
772  Fonction de fab pour creer le slaving.
773  La fonction cree un tableau de booleens avec une seule spider.
774  Utilisee par la fonction compute6Segments()
775  """
776  a = np.ones((N, N), dtype=bool)
777  X = (np.arange(N) - i0) * scale
778  Y = (np.arange(N) - j0) * scale
779  X, Y = np.meshgrid(X, Y, indexing='ij') # convention d'appel [x,y]
780  w = 2 * np.pi / 6
781  i = nspider
782  nn = (abs(X * np.cos(i * w - rot) + Y * np.sin(i * w - rot)) < dspider / 2.)
783  a[nn] = False
784  return a
785 
786 
787 def fillSpider(N, nspider, dspider, i0, j0, scale, rot):
788  """
789  Creates a boolean spider mask on a map of dimensions (N,N)
790  The spider is centred at floating-point coords (i0,j0).
791 
792  :returns: spider image (boolean)
793  :param int N: size of output image
794  :param int nspider: number of spiders
795  :param float dspider: width of spiders
796  :param float i0: coord of spiders symmetry centre
797  :param float j0: coord of spiders symmetry centre
798  :param float scale: size of a pixel in same unit as dspider
799  :param float rot: rotation angle in radians
800 
801  """
802  a = np.ones((N, N), dtype=bool)
803  X = (np.arange(N) - i0) * scale
804  Y = (np.arange(N) - j0) * scale
805  X, Y = np.meshgrid(X, Y, indexing='ij') # convention d'appel [x,y]
806  w = 2 * np.pi / nspider
807  for i in range(nspider):
808  nn = (abs(X * np.cos(i * w - rot) + Y * np.sin(i * w - rot)) < dspider / 2.)
809  a[nn] = False
810  return a
811 
812 
813 def fillHalfSpider(N, nspider, dspider, i0, j0, scale, rot):
814  a = np.ones((N, N), dtype=bool)
815  b = np.ones((N, N), dtype=bool)
816  X = (np.arange(N) - i0) * scale
817  Y = (np.arange(N) - j0) * scale
818  X, Y = np.meshgrid(X, Y, indexing='ij') # convention d'appel [x,y]
819  w = 2 * np.pi / nspider
820  for i in range(nspider):
821  right = (X * np.cos(i * w - rot) + Y * np.sin(i * w - rot) <
822  dspider / 2) * (X * np.cos(i * w - rot) + Y * np.sin(i * w - rot) > 0.)
823  left = (X * np.cos(i * w - rot) + Y * np.sin(i * w - rot) >
824  -dspider / 2) * (X * np.cos(i * w - rot) + Y * np.sin(i * w - rot) < 0.)
825  a[right] = False
826  b[left] = False
827  return a, b
828 
829 
830 def createHexaPattern(pitch, supportSize):
831  """
832  Cree une liste de coordonnees qui decrit un maillage hexagonal.
833  Retourne un tuple (x,y).
834 
835  Le maillage est centre sur 0, l'un des points est (0,0).
836  Une des pointes de l'hexagone est dirigee selon l'axe Y, au sens ou le
837  tuple de sortie est (x,y).
838 
839  :param float pitch: distance between 2 neighbour points
840  :param int supportSize: size of the support that need to be populated
841 
842  """
843  V3 = np.sqrt(3)
844  nx = int(np.ceil((supportSize / 2.0) / pitch) + 1)
845  x = pitch * (np.arange(2 * nx + 1) - nx)
846  ny = int(np.ceil((supportSize / 2.0) / pitch / V3) + 1)
847  y = (V3 * pitch) * (np.arange(2 * ny + 1) - ny)
848  x, y = np.meshgrid(x, y, indexing='ij')
849  x = x.flatten()
850  y = y.flatten()
851  peak_axis = np.append(x, x + pitch / 2.) # axe dirige selon sommet
852  flat_axis = np.append(y, y + pitch * V3 / 2.) # axe dirige selon plat
853  return flat_axis, peak_axis
854 
855 
856 def generateCoordSegments(D, rot, pitch=1.244683637214, nseg=33, inner_rad=4.1,
857  outer_rad=15.4, R=95.7853, nominalD=40):
858  """
859  Computes the coordinates of the corners of all the hexagonal
860  segments of M1.
861  Result is a tuple of arrays(6, 798).
862 
863  Args:
864  D: (float) : pupil diameter in meters (it must be set to 40.0 m for the ELT)
865 
866  rot: (float) : pupil rotation angle in radians
867 
868  pitch: (float): Segment pitch [meters]
869 
870  nseg: (int) : number of segments across the diameter
871 
872  inner_rad : (float): Inner radius [meters]
873 
874  outer_rad : (float): Outer radius [meters]
875 
876  R : (float): Curvature radius of the M1
877 
878  nominalD: (float): diameter for nominal pupil
879 
880  """
881  V3 = np.sqrt(3)
882  #pitch = 1.227314 # no correction du bol
883  #pitch = 1.244683637214 # diametre du cerle INSCRIT
884  # diamseg = pitch*2/V3 # diametre du cercle contenant TOUT le segment
885  # print("segment diameter : %.6f\n" % diamseg)
886 
887  # Creation d'un pattern hexa avec pointes selon la variable <ly>
888  lx, ly = createHexaPattern(pitch, (nseg + 2) * pitch)
889  ll = np.sqrt(lx**2 + ly**2)
890  # Elimination des segments non valides grace a 2 nombres parfaitement
891  # empiriques ajustes a-la-mano.
892  #inner_rad, outer_rad = 4.1, 15.4 # nominal, 798 segments
893  nn = (ll > inner_rad * pitch) & (ll < outer_rad * pitch)
894  lx = lx[nn]
895  ly = ly[nn]
896  lx, ly = reorganizeSegmentsOrderESO(lx, ly)
897  ll = np.sqrt(lx**2 + ly**2)
898 
899  # n = ll.shape[0]
900  # print("Nbre de segments : %d\n" % n)
901  # Creation d'un hexagone-segment avec pointe dirigee vers
902  # variable <hx> (d'ou le cos() sur hx)
903  th = np.linspace(0, 2 * np.pi, 7)[0:6]
904  hx = np.cos(th) * pitch / V3
905  hy = np.sin(th) * pitch / V3
906 
907  # Le maillage qui permet d'empiler des hexagones avec sommets 3h-9h
908  # est un maillage hexagonal avec sommets 12h-6h, donc a 90°.
909  # C'est pour ca qu'il a fallu croiser les choses avant.
910  x = (lx[None, :] + hx[:, None])
911  y = (ly[None, :] + hy[:, None])
912  r = np.sqrt(x**2 + y**2)
913  #R = 95.7853
914  rrc = R / r * np.arctan(r / R) # correction factor
915  x *= rrc
916  y *= rrc
917 
918  #nominalD = 40.0 # size of the OFFICIAL E-ELT
919  if D != nominalD:
920  x *= D / nominalD
921  y *= D / nominalD
922 
923  # Rotation matrices
924  mrot = np.array([[np.cos(rot), np.sin(rot)], [-np.sin(rot), np.cos(rot)]])
925 
926  # rotation of coordinates
927  # le tableau [x,y] est de taille (2,6,798). Faut un transpose a la con
928  # pour le transformer en (6,2,798) pour pouvoir faire le np.dot
929  # correctement. En sortie, xrot est (2,6,798).
930  xyrot = np.dot(mrot, np.transpose(np.array([x, y]), (1, 0, 2)))
931 
932  return xyrot[0], xyrot[1]
933 
934 
935 def gendron():
936  """
937  La fonction est appelee quand l'utilisateur a demande une pupille
938  ELT, et renseigne un diametre de telescope different de 40 metres.
939 
940  Faut vraiment que je commente ou t'as compris ??
941 
942  """
943  mymsg = [
944  "\n\n\n\n",
945  # "__ ___ ____ _ _ ___ _ _ ___ _",
946  # "\ \ / / \ | _ \| \ | |_ _| \ | |/ ___|",
947  # " \ \ /\ / / _ \ | |_) | \| || || \| | | _ ",
948  # " \ V V / ___ \| _ <| |\ || || |\ | |_| |",
949  # " \_/\_/_/ \_\_| \_\_| \_|___|_| \_|\____|", " \n",
950  "Vous utilisez un telescope de type ELT. Ce telescope",
951  "est fait pour etre utilise avec un diametre de 40 m.", " ",
952  "Or, vous utilisez un diametre different. Cela signifie",
953  "que le telescope que vous etes en train de creer a une",
954  "taille differente du veritable E-ELT de l'ESO.", " ",
955  " * Soit vous savez exactement ce que vous faites, auquel",
956  "cas bonne route.", " ",
957  " * Soit vous desirez creer LE vrai E-ELT et il faut changer",
958  "plusieurs choses:",
959  " 1) le diametre telescope de votre fichier de parametres et",
960  " le renseigner a 40 metres.",
961  " p_tel.set_diam(40.0) # Nominal size for the real EELT",
962  " 2) le nombre d'actionneurs de M4 a 75",
963  " p_dm0.set_nact(75) # 75 actu in 40m for pitch=54.05cm",
964  " 3) option: tourner la pupille de 90 degres pour revenir au",
965  " cas initial de compass",
966  " p_tel.set_pupangle(90.) # ELT pup rotation in degrees"
967  " ", "\n\n"
968  ]
969  for ligne in mymsg:
970  print(ligne)
971 
972 
973 def reorganizeSegmentsOrderESO(x, y):
974  """
975  Reorganisation des segments facon ESO.
976  Voir
977  ESO-193058 Standard Coordinate System and Basic Conventions
978 
979  :param float x: tableau des centres X des segments
980  :param float y: idem Y
981  :return tuple (x,y): meme tuple que les arguments d'entree, mais tries.
982 
983  """
984  # pi/2, pi/6, 2.pi, ...
985  pi_3 = np.pi / 3
986  pi_6 = np.pi / 6
987  pix2 = 2 * np.pi
988  # calcul des angles
989  t = (np.arctan2(y, x) + pi_6 - 1e-3) % (pix2)
990  X = np.array([])
991  Y = np.array([])
992  A = 100.
993  for k in range(6):
994  sector = (t > k * pi_3) & (t < (k + 1) * pi_3)
995  u = k * pi_3
996  distance = (A * np.cos(u) - np.sin(u)) * x[sector] + (
997  np.cos(u) + A * np.sin(u)) * y[sector]
998  indsort = np.argsort(distance)
999  X = np.append(X, x[sector][indsort])
1000  Y = np.append(Y, y[sector][indsort])
1001  return X, Y
1002 
1003 
1004 def getdatatype(truc):
1005  """
1006  Returns the data type of a numpy variable, either scalar value or array
1007  """
1008  if np.isscalar(truc):
1009  return type(truc)
1010  else:
1011  return type(truc.flatten()[0])
1012 
1013 
1014 def generateSegmentProperties(attribute, hx, hy, i0, j0, scale, gap, N, D, softGap=0,
1015  nominalD=40, pitch=1.244683637214, half_seg=0.75):
1016  """
1017  Builds a 2D image of the pupil with some attributes for each of the
1018  segments. Those segments are described from arguments hx and hy, that
1019  are produced by the function generateCoordSegments(D, rot).
1020 
1021  When attribute is a phase, then it must be a float array of dimension
1022  [3, 798] with the dimension 3 being piston, tip, and tilt.
1023  Units of phase is xxx rms, and the output of the procedure will be
1024  in units of xxx.
1025 
1026 
1027  :returns: pupil image (N, N), with the same type of input argument attribute
1028 
1029  :param float/int/bool attribute: scalar value or 1D-array of the reflectivity of
1030  the segments or 2D array of phase
1031  If attribute is scalar, the value will be replicated for all segments.
1032  If attribute is a 1D array, then it shall contain the reflectivities
1033  of all segments.
1034  If attribute is a 2D array then it shall contain the piston, tip
1035  and tilt of the segments. The array shall be of dimension
1036  [3, 798] that contains [piston, tip, tilt]
1037  On output, the data type of the pupil map will be the same as attribute.
1038  :param float hx, hy: arrays [6,:] describing the segment shapes. They are
1039  generated using generateCoordSegments()
1040  :param float dspider: width of spiders in meters
1041  :param float i0, j0: index of pixels where the pupil should be centred.
1042  Can be floating-point indexes.
1043  :param float scale: size of a pixel of the image, in meters.
1044  :param float gap: half-space between segments in meters
1045  :param int N: size of the output array (N,N)
1046  :param float D: diameter of the pupil. For the nominal EELT, D shall
1047  be set to 40.0
1048  :param bool softGap: if False, the gap between segments is binary 0/1
1049  depending if the pixel is within the gap or not. If True, the gap
1050  is a smooth region of a fwhm of 2 pixels with a depth related to the
1051  gap width.
1052  :param float nominalD: diameter needed to get nominal pupil aperture
1053  :param float pitch: segment pitch
1054  :param float half_seg: segment half size
1055 
1056 
1057 
1058  attribute = np.ones(798)+np.random.randn(798)/20.
1059  N = 800
1060  i0 = N/2
1061  j0 = N/2
1062  rotdegree = 0.0
1063  scale = 41./N
1064  gap = 0.03
1065 
1066  """
1067 
1068  # number of segments
1069  nseg = hx.shape[-1]
1070  # If <attribute> is a scalar, then we make a list. It will be required
1071  # later on to set the attribute to each segment.
1072  if np.isscalar(attribute):
1073  attribute = np.array([attribute] * nseg)
1074 
1075  # the pupil map is created with the same data type as <attribute>
1076  pupil = np.zeros((N, N), dtype=getdatatype(attribute))
1077 
1078  # average coord of segments
1079  x0 = np.mean(hx, axis=0)
1080  y0 = np.mean(hy, axis=0)
1081  # avg coord of segments in pixel indexes
1082  x0 = x0 / scale + i0
1083  y0 = y0 / scale + j0
1084  # size of mini-support
1085  hexrad = half_seg * D / nominalD / scale
1086  ix0 = np.floor(x0 - hexrad).astype(int) - 1
1087  iy0 = np.floor(y0 - hexrad).astype(int) - 1
1088  segdiam = np.ceil(hexrad * 2 + 1).astype(int) + 1
1089 
1090  n = attribute.shape[0]
1091  if n != 3:
1092  # attribute is a signel value : either reflectivity, or boolean,
1093  # or just piston.
1094  if softGap != 0:
1095  # Soft gaps
1096  # The impact of gaps are modelled using a simple function: Lorentz, 1/(1+x**2)
1097  # The fwhm is always equal to 2 pixels because the gap is supposed
1098  # to be "small/invisible/undersampled". The only visible thing is
1099  # the width of the impulse response, chosen 2-pixel wide to be
1100  # well sampled.
1101  # The "depth" is related to the gap width. The integral of a Lorentzian
1102  # of 2 pix wide is PI. Integral of a gap of width 'gap' in pixels is 'gap'.
1103  # So the depth equals to gap/scale/np.pi.
1104  for i in range(nseg):
1105  indx, indy, distedge = fillPolygon(hx[:, i], hy[:, i], i0 - ix0[i],
1106  j0 - iy0[i], scale, gap * 0., segdiam,
1107  index=1)
1108  pupil[indx + ix0[i], indy + iy0[i]] = attribute[i] * (
1109  1. - (gap / scale / np.pi) / (1 + (distedge / scale)**2))
1110  else:
1111  # Hard gaps
1112  for i in range(nseg):
1113  indx, indy, distedge = fillPolygon(hx[:, i], hy[:, i], i0 - ix0[i],
1114  j0 - iy0[i], scale, gap, segdiam,
1115  index=1)
1116  pupil[indx + ix0[i], indy + iy0[i]] = attribute[i]
1117  else:
1118  # attribute is [piston, tip, tilt]
1119  minimap = np.zeros((segdiam, segdiam))
1120  xmap = np.arange(segdiam) - segdiam / 2
1121  xmap, ymap = np.meshgrid(xmap, xmap, indexing='ij') # [x,y] convention
1122  #pitch = 1.244683637214 # diameter of inscribed circle
1123  diamseg = pitch * 2 / np.sqrt(3) # diameter of circumscribed circle
1124  diamfrizou = (pitch + diamseg) / 2 * D / nominalD # average diameter of the 2
1125  # Calcul du facteur de mise a l'echelle pour l'unite des tilts.
1126  # xmap et ymap sont calculees avec un increment de +1 pour deux pixels
1127  # voisins, donc le facteur a appliquer est tel que l'angle se conserve
1128  # donc factunit*1 / scale = 4*factunit
1129  factunit = 4 * scale / diamfrizou
1130  for i in range(nseg):
1131  indx, indy, _ = fillPolygon(hx[:, i], hy[:, i], i0 - ix0[i], j0 - iy0[i],
1132  scale, 0., segdiam, index=1)
1133  minimap = attribute[0, i] + (factunit * attribute[1, i]) * xmap + (
1134  factunit * attribute[2, i]) * ymap
1135  pupil[indx + ix0[i], indy + iy0[i]] = minimap[indx, indy]
1136 
1137  return pupil
1138 
1139 
1140 def centrePourVidal(N, i0, j0, centerMark):
1141  """
1142  Renvoie une image de boolens (False) de taille (N,N) avec un point
1143  ou une croix (True) centree sur (i0, j0).
1144  :param int N: taille de l'image de sortie
1145  :param float i0, j0: position du marqueur de sortie
1146  :param int centerMark: 0 (pour rien), 1 (option point) ou 2 (option croix)
1147  """
1148  scale = 1.0
1149  res = 0
1150  X = (np.arange(N) - i0) * scale
1151  Y = (np.arange(N) - j0) * scale
1152  X, Y = np.meshgrid(X, Y, indexing='ij') # convention d'appel [x,y]
1153  if centerMark == 1:
1154  res = (X**2 + Y**2) < 1
1155  if centerMark == 2:
1156  res = (np.abs(X) < 0.9) | (np.abs(Y) < 0.9)
1157  return res
shesha.constants
Numerical constants for shesha and config enumerations for safe-typing.
Definition: constants.py:1
shesha.util.make_pupil.fillPolygon
def fillPolygon(x, y, i0, j0, scale, gap, N, index=0)
From a list of points defined by their 2 coordinates list x and y, creates a filled polygon with side...
Definition: make_pupil.py:661
shesha.util.make_pupil.centrePourVidal
def centrePourVidal(N, i0, j0, centerMark)
Renvoie une image de boolens (False) de taille (N,N) avec un point ou une croix (True) centree sur (i...
Definition: make_pupil.py:1155
shesha.util.make_pupil.generateEeltPupilMask
def generateEeltPupilMask(npt, dspider, i0, j0, pixscale, gap, rotdegree, D=40.0, cobs=0, centerMark=0, halfSpider=False, pitch=1.244683637214, nseg=33, inner_rad=4.1, outer_rad=15.4, R=95.7853, nominalD=40, half_seg=0.75, refl=None, rotSpiderDegree=None, nmissing=0)
Generates a boolean pupil mask of the binary EELT pupil on a map of size (npt, npt).
Definition: make_pupil.py:555
shesha.util.make_pupil.reorganizeSegmentsOrderESO
def reorganizeSegmentsOrderESO(x, y)
Reorganisation des segments facon ESO.
Definition: make_pupil.py:989
shesha.util.make_pupil.make_VLT
def make_VLT(dim, pupd, tel)
Initialize the VLT pupil.
Definition: make_pupil.py:218
shesha.util.make_pupil.make_phase_ab
def make_phase_ab(dim, pupd, tel, pup=None, xc=-1, yc=-1, real=0, halfSpider=False)
Compute the EELT M1 phase aberration.
Definition: make_pupil.py:386
shesha.util.make_pupil.compute6Segments
def compute6Segments(pupNoSpiders, N, pixscale, dspider, i0, j0, rot=0)
Definition: make_pupil.py:748
shesha.util.make_pupil.createHexaPattern
def createHexaPattern(pitch, supportSize)
Cree une liste de coordonnees qui decrit un maillage hexagonal.
Definition: make_pupil.py:848
shesha.util.make_pupil.make_pupil
def make_pupil(dim, pupd, tel, xc=-1, yc=-1, real=0, halfSpider=False)
Initialize the system pupil.
Definition: make_pupil.py:68
shesha.util.make_pupil.make_EELT
def make_EELT(dim, pupd, tel, N_seg=-1)
Initialize the EELT pupil.
Definition: make_pupil.py:265
shesha.util.make_pupil.generateCoordSegments
def generateCoordSegments(D, rot, pitch=1.244683637214, nseg=33, inner_rad=4.1, outer_rad=15.4, R=95.7853, nominalD=40)
Computes the coordinates of the corners of all the hexagonal segments of M1.
Definition: make_pupil.py:886
shesha.util.make_pupil.gendron
def gendron()
La fonction est appelee quand l'utilisateur a demande une pupille ELT, et renseigne un diametre de te...
Definition: make_pupil.py:948
shesha.util.make_pupil.compute1Spider
def compute1Spider(nspider, N, dspider, i0, j0, scale, rot)
Fonction de fab pour creer le slaving.
Definition: make_pupil.py:781
shesha.util.make_pupil.getdatatype
def getdatatype(truc)
Returns the data type of a numpy variable, either scalar value or array.
Definition: make_pupil.py:1013
shesha.util.make_pupil.fillHalfSpider
def fillHalfSpider(N, nspider, dspider, i0, j0, scale, rot)
Definition: make_pupil.py:819
shesha.util.make_pupil.generateSegmentProperties
def generateSegmentProperties(attribute, hx, hy, i0, j0, scale, gap, N, D, softGap=0, nominalD=40, pitch=1.244683637214, half_seg=0.75)
Builds a 2D image of the pupil with some attributes for each of the segments.
Definition: make_pupil.py:1074
shesha.util.make_pupil.fillSpider
def fillSpider(N, nspider, dspider, i0, j0, scale, rot)
Creates a boolean spider mask on a map of dimensions (N,N) The spider is centred at floating-point co...
Definition: make_pupil.py:807
shesha.util.make_pupil.make_pupil_generic
def make_pupil_generic(dim, pupd, t_spiders=0.01, spiders_type=SpiderType.SIX, xc=0, yc=0, real=0, cobs=0)
Initialize the system pupil.
Definition: make_pupil.py:148