39 from astropy.io
import fits
44 from .
import utilities
as util
46 from typing
import List, Union
50 """ Compute the DM support dimensions
54 cent : (float): center of the pupil
56 extent: (float): size of the DM support
58 ssize: (int): size of ipupil support
60 n1 = np.floor(cent - extent / 2)
61 n2 = np.ceil(cent + extent / 2)
67 return int(n1), int(n2)
70 def dim_dm_patch(pupdiam: int, diam: float, type: bytes, alt: float,
71 xpos_wfs: List[float], ypos_wfs: List[float]):
72 """ compute patchDiam for DM
76 pupdiam: (int) : pupil diameter
78 diam: (float) : telescope diameter
80 type: (bytes) : type of dm
82 alt: (float) : altitude of dm
84 xpos_wfs: (list) : list of wfs xpos
86 ypos_wfs: (list) : list of wfs ypos
89 if len(xpos_wfs) == 0:
93 np.linalg.norm([xpos_wfs[w], ypos_wfs[w]])
for w
in range(len(xpos_wfs))
95 if ((type == scons.DmType.PZT)
or (type == scons.DmType.TT)):
97 elif (type == scons.DmType.KL):
100 raise TypeError(
"This type of DM doesn't exist ")
102 patchDiam = int(pupdiam + 2 * np.max(norms) * CONST.ARCSEC2RAD * np.abs(alt) / (pp))
108 Creates a list of M=nxact^2 actuator positions spread over an square grid.
109 Coordinates are centred around (0,0).
113 pitch: (float) : distance in pixels between 2 adjacent actus
115 nxact: (int) : number of actu across the pupil diameter
119 xy: (np.ndarray(dims=2,dtype=np.float32)) : xy[M,2] list of coodinates
122 xy = np.tile(np.arange(nxact) - (nxact - 1.) / 2., (nxact, 1)).astype(np.float32)
123 xy = np.array([xy.flatten(), xy.T.flatten()]) * pitch
130 Creates a list of M actuator positions spread over an hexagonal grid.
131 The number M is the number of points of this grid, it cannot be
132 known before the procedure is called.
133 Coordinates are centred around (0,0).
134 The support that limits the grid is a square [-supportSize/2, supportSize/2].
138 pitch: (float) : distance in pixels between 2 adjacent actus
140 supportSize: (int) : size in pixels of the support over which the coordinate list
145 xy: (np.ndarray(dims=2,dtype=np.float32)) : xy[2,M] list of coordinates
148 nx = int(np.ceil((supportSize / 2.0) / pitch) + 1)
149 x = pitch * (np.arange(2 * nx + 1, dtype=np.float32) - nx)
151 ny = int(np.ceil((supportSize / 2.0) / pitch / V3) + 1)
152 y = (V3 * pitch) * (np.arange(2 * ny + 1, dtype=np.float32) - ny)
154 x = np.tile(x, (Ny, 1)).flatten()
155 y = np.tile(y, (Nx, 1)).T.flatten()
156 x = np.append(x, x + pitch / 2.)
157 y = np.append(y, y + pitch * V3 / 2.)
158 xy = np.float32(np.array([y, x]))
164 Creates a list of M actuator positions spread over an hexagonal grid.
165 The number M is the number of points of this grid, it cannot be
166 known before the procedure is called.
167 Coordinates are centred around (0,0).
168 The support of the grid is a square [-supportSize/2,vsupportSize/2].
172 pitch: (float) : distance in pixels between 2 adjacent actus
173 supportSize: (int) : size in pixels of the support over which the coordinate list
175 pupAngleDegree: (float) : Rotation angle of the DM
179 xy: (np.ndarray(dims=2,dtype=np.float32)) : xy[2,M] list of coodinates
185 ny = int(np.ceil((supportSize / 2.0) / pitch) + 1)
186 y = pitch * (np.arange(2 * ny + 1, dtype=np.float32) - ny)
187 nx = int(np.ceil((supportSize / 2.0) / pitch / V3) + 1)
188 x = (V3 * pitch) * (np.arange(2 * nx + 1, dtype=np.float32) - nx)
191 x, y = np.meshgrid(x, y, indexing=
'ij')
194 x = np.append(x, x + pitch * V3 / 2.)
195 y = np.append(y, y + pitch / 2.)
204 th = np.arctan2(y, x)
205 nn = np.where(((th < pi / 6) & (th > -pi / 6)))
214 xx = np.cos(k * pi / 3) * x + np.sin(k * pi / 3) * y
215 yy = -np.sin(k * pi / 3) * x + np.cos(k * pi / 3) * y
220 rot = pupAngleDegree * np.pi / 180.0
221 mrot = np.array([[np.cos(rot), -np.sin(rot)], [np.sin(rot), np.cos(rot)]])
222 XY = np.dot(mrot, [X, Y])
223 return np.float32(XY)
227 threshold: float) -> np.ndarray:
229 Select actuators based on their distance to the nearest pupil pixel
230 The implementation proposed here is easy but limits it precision to
231 an integer roundoff of the threshold
233 actuPos: 2 x nActu np.array[float]: actuator position list - pupil pixel units
234 pupil: nPup x nPup np.ndarray[bool]: pupil mask
235 threshold: float: max allowed distance - pupil pixel units
239 from scipy.ndimage.morphology
import binary_dilation
240 k = np.ceil(threshold)
241 i, j = np.meshgrid(np.arange(2 * k + 1), np.arange(2 * k + 1), indexing=
'ij')
242 disk = ((i - k)**2 + (j - k)**2)**.5 <= k
243 dilatedPupil = binary_dilation(pupil, disk)
244 actuIsIn = dilatedPupil[(np.round(actuPos[0]).astype(np.int32),
245 np.round(actuPos[1]).astype(np.int32))]
246 return actuPos[:, actuIsIn]
249 def select_actuators(xc: np.ndarray, yc: np.ndarray, nxact: int, pitch: int, cobs: float,
250 margin_in: float, margin_out: float, N=
None):
252 Select the "valid" actuators according to the system geometry
256 xc: actuators x positions (origine in center of mirror)
258 yc: actuators y positions (origine in center of mirror)
274 liste_fin: actuator indice selection for xpos/ypos
280 dis = np.sqrt(xc**2 + yc**2)
283 rad_in = (((nxact - 1) / 2) * cobs - margin_in) * pitch
286 if (margin_out
is None):
288 rad_out = ((nxact - 1.) / 2. + margin_out) * pitch
290 valid_actus = np.where((dis <= rad_out) * (dis >= rad_in))[0]
293 valid_actus = np.where(dis >= rad_in)[0]
294 indsort = np.argsort(dis[valid_actus])
296 if (N > valid_actus.size):
297 print(
'Too many actuators wanted, restricted to ', valid_actus.size)
299 valid_actus = np.sort(indsort[:N])
304 def make_zernike(nzer: int, size: int, diameter: int, xc=-1., yc=-1., ext=0):
305 """Compute the zernike modes
309 nzer: (int) : number of modes
311 size: (int) : size of the screen
313 diameter: (int) : pupil diameter
315 xc: (float) : (optional) x-position of the center
317 yc: (float) : (optional) y-position of the center
319 ext: (int) : (optional) extension
323 z : (np.ndarray(ndims=3,dtype=np.float64)) : zernikes modes
333 radius = (diameter) / 2.
334 zr = util.dist(size, xc, yc).astype(np.float32).T / radius
335 zmask = np.zeros((zr.shape[0], zr.shape[1], nzer), dtype=np.float32)
336 zmaskmod = np.zeros((zr.shape[0], zr.shape[1], nzer), dtype=np.float32)
338 zmask[:, :, 0] = (zr <= 1).astype(np.float32)
339 zmaskmod[:, :, 0] = (zr <= 1.2).astype(np.float32)
341 for i
in range(1, nzer):
342 zmask[:, :, i] = zmask[:, :, 0]
343 zmaskmod[:, :, i] = zmaskmod[:, :, 0]
345 zrmod = zr * zmaskmod[:, :, 0]
347 zr = zr * zmask[:, :, 0]
349 x = np.tile(np.arange(size).astype(np.float32), (size, 1))
350 zteta = np.arctan2(x - yc, x.T - xc).astype(np.float32)
352 z = np.zeros((size, size, nzer), dtype=np.float32)
354 for zn
in range(nzer):
358 for i
in range((n - m) // 2 + 1):
359 z[:, :, zn] = z[:, :, zn] + (-1.) ** i * zrmod ** (n - 2. * i) * float(np.math.factorial(n - i)) / \
360 float(np.math.factorial(i) * np.math.factorial((n + m) // 2 - i) *
361 np.math.factorial((n - m) // 2 - i))
363 for i
in range((n - m) // 2 + 1):
364 z[:, :, zn] = z[:, :, zn] + (-1.) ** i * zr ** (n - 2. * i) * float(np.math.factorial(n - i)) / \
365 float(np.math.factorial(i) * np.math.factorial((n + m) // 2 - i) *
366 np.math.factorial((n - m) // 2 - i))
368 if ((zn + 1) % 2 == 1):
370 z[:, :, zn] = z[:, :, zn] * np.sqrt(n + 1.)
372 z[:, :, zn] = z[:, :, zn] * \
373 np.sqrt(2. * (n + 1)) * np.sin(m * zteta)
376 z[:, :, zn] = z[:, :, zn] * np.sqrt(n + 1.)
378 z[:, :, zn] = z[:, :, zn] * \
379 np.sqrt(2. * (n + 1)) * np.cos(m * zteta)
389 Returns the radial degree and the azimuthal number of zernike
390 number zn, according to Noll numbering (Noll, JOSA, 1976)
394 zn: (int) : zernike number
398 rd: (int) : radial degrees
400 an: (int) : azimuthal numbers
405 for m
in range(n + 1):
406 if ((n - m) % 2 == 0):
416 dm_fits_content=
"""The DM FITS file is compatible with COMPASS DM database.
417 The primary header contains the keywords:
418 * PIXSIZE : the size of the pixels on the maps in meters.
420 * XCENTER, YCENTER are the coordinates of the centre of the pupil, expressed in pixels, in a reference frame conformable to (i,j) coords. The translation from pixels to meters can be done using:
421 meters = (pixels - XCENTER) * PIXSIZE
423 * PUPM is the diameter of pupil stop (meters).
425 * Additionally the header provides the user with:
426 PITCHM is the size of the DM pitch in meters (may be handy is some cases and useful when using this file with COMPASS software)
428 This FITS file contains 3 extensions:
429 * Extension 'I1_J1' are the coordinate (i,j) of the first pixel for each of the 2D maps (see Extension 2), so that they can be inserted in a larger map.
431 * Extension 'INFLU' are the 2D maps of the influence functions.
433 * Extension 'XPOS_YPOS' are the coordinates (xpos, ypos) of the physical location of the actuator, in pixels. This data is provided for information only and does not directly participate to build the DM. The present coordinates are positions in M1 space, i.e. include the distorsion due to telescope optics.
437 """adds content to a docstring (to be used as decorator)"""
439 obj.__doc__ = obj.__doc__.format(content)
443 def write_dm_custom_fits(file_name, i1, j1, influ_cube, xpos, ypos, xcenter, ycenter,pixsize,pupm, *,pitchm=None):
444 """Write a custom_dm fits file based on user provided data (see args)
447 file_name : (string) : name of the custom dm fits file
449 i1 : (np.ndarray) : x coordinate of the first pixel for each of the 2D maps
451 j1 : (np.ndarray) : y coordinate of the first pixel for each of the 2D maps
453 influ_cube : (np.ndarray) : 2D maps of the influence functions
455 xpos : (np.ndarray) : x coordinate of the physical location of the actuator
457 ypos : (np.ndarray) : y coordinate of the physical location of the actuator
459 xcenter : (float) : x coordinate of the centre of the pupil, expressed in pixels
461 ycenter : (float) : y coordinate of the centre of the pupil, expressed in pixels
463 pixsize : (float) : size of the pixels on the maps in meters
465 pupm : (float) : diameter of pupil stop (meters)
468 pitchm : (float) : size of the DM pitch in meters. Defaults to None.
471 (HDUList) : custom_dm data
474 primary_hdu = fits.PrimaryHDU()
475 primary_hdu.header[
'VERSION'] = (fits_version,
'file format version')
476 primary_hdu.header[
'XCENTER'] = (xcenter ,
'DM centre along X in pixels')
477 primary_hdu.header[
'YCENTER'] = (ycenter ,
'DM centre along Y in pixels')
478 primary_hdu.header[
'PIXSIZE'] = (pixsize ,
'pixel size (meters)')
479 primary_hdu.header[
'PUPM'] = (pupm ,
'nominal pupil diameter (meters)')
480 if(pitchm
is not None):
481 primary_hdu.header[
'PITCHM'] = (pitchm,
'DM pitch (meters)')
483 for line
in dm_fits_content.splitlines():
484 primary_hdu.header.add_comment(line)
486 image_hdu = fits.ImageHDU(np.c_[i1 , j1 ].T, name=
"I1_J1")
487 image_hdu2 = fits.ImageHDU(influ_cube, name=
"INFLU")
488 image_hdu3 = fits.ImageHDU(np.c_[xpos, ypos].T, name=
"XPOS_YPOS")
490 dm_custom = fits.HDUList([primary_hdu, image_hdu, image_hdu2, image_hdu3])
492 dm_custom.writeto(file_name,overwrite=1)
495 @add_doc_content(dm_fits_content)
497 """Return an HDUList (FITS) with the data required to create a COMPASS custom_dm
502 p_dm : (Param_dm) : dm settings
504 p_geom : (Param_geom) : geometry settings
507 file_name : (string) : if set, the HDU is written to the file specified by this variable
509 p_tel : (Param_tel) : telescope settings, used to provide the diameter (if not provided, the default diameter id obtained from the p_geom as pupdiam*pixsize)
512 (HDUList) : custom_dm data
515 pixsize = p_geom.get_pixsize()
516 diam = p_geom.pupdiam * p_geom._pixsize
517 if(p_tel
is not None):
521 i1 = p_dm._i1 + p_dm._n1
522 j1 = p_dm._j1 + p_dm._n1
523 influ = p_dm._influ / p_dm.unitpervolt
525 xcenter = p_geom.cent
526 ycenter = p_geom.cent
530 pitchm = p_dm._pitch*pixsize
531 dm_custom =
write_dm_custom_fits(file_name,i1,j1,influ,xpos,ypos,xcenter,ycenter,pixsize,diam,pitchm=pitchm)
Numerical constants for shesha and config enumerations for safe-typing.
def createDoubleHexaPattern(float pitch, int supportSize, float pupAngleDegree)
Creates a list of M actuator positions spread over an hexagonal grid.
def dim_dm_support(float cent, int extent, int ssize)
Compute the DM support dimensions.
def select_actuators(np.ndarray xc, np.ndarray yc, int nxact, int pitch, float cobs, float margin_in, float margin_out, N=None)
Select the "valid" actuators according to the system geometry.
def export_custom_dm(file_name, p_dm, p_geom, *p_tel=None)
Return an HDUList (FITS) with the data required to create a COMPASS custom_dm.
def write_dm_custom_fits(file_name, i1, j1, influ_cube, xpos, ypos, xcenter, ycenter, pixsize, pupm, *pitchm=None)
Write a custom_dm fits file based on user provided data (see args)
def zernumero(int zn)
Returns the radial degree and the azimuthal number of zernike number zn, according to Noll numbering ...
def make_zernike(int nzer, int size, int diameter, xc=-1., yc=-1., ext=0)
Compute the zernike modes.
def add_doc_content(*content)
adds content to a docstring (to be used as decorator)
np.ndarray filterActuWithPupil(np.ndarray actuPos, np.ndarray pupil, float threshold)
Select actuators based on their distance to the nearest pupil pixel The implementation proposed here ...
def dim_dm_patch(int pupdiam, float diam, bytes type, float alt, List[float] xpos_wfs, List[float] ypos_wfs)
compute patchDiam for DM
def createHexaPattern(float pitch, int supportSize)
Creates a list of M actuator positions spread over an hexagonal grid.
def createSquarePattern(float pitch, int nxact)
Creates a list of M=nxact^2 actuator positions spread over an square grid.