Optical functions

Propagation functions

Asterix.optics.propagation_functions.mft(image, real_dim_input=4, dim_output=4, nbres=1, inverse=False, norm='backward', X_offset_input=0, Y_offset_input=0, X_offset_output=0, Y_offset_output=0, only_mat_mult=False, AA=None, BB=None, norm0=None, returnAABB=False, dtype_complex='complex128')

Return the Matrix Direct Fourier transform (MFT) of a 2D image.

Based on Matrix Direct Fourier transform (MFT) from R. Galicher (cf.Soummer et al. 2007, OSA). Can deal with any size, any position of the 0-frequency…

This function measures 2 matrices AA and BB, and the normalization factor ‘norm0’. The MFT itself is the matrix multiplication norm0 * ((AA @ image) @ BB) (where @ is a matrix multiplication)

Can be used in a classical way: ZeMFT = mft(image,

real_dim_input=real_dim_input, dim_output=dim_output, nbres=nbres, only_mat_mult=False, returnAABB=False)

or you can measure AA, BB, and norm0 and / or only do the matrix multiplication: AA,BB,norm0 = mft(image,

real_dim_input=real_dim_input, dim_output=dim_output, nbres=nbres, only_mat_mult=False, returnAABB=True)

ZeMFT = mft(image,

AA=AA, BB=BB, norm0 = norm0, only_mat_mult=True, returnAABB=False)

By separating those 2 steps you can save a lot of time. If you are doing a lot of MFTs with the same input and output dimension parameters, only the second step (with only_mat_mult=True) need to be done.

AUTHORS: Baudoz, Galicher, Mazoyer

REVISION HISTORY :

-Revision 1.1 2011 Initial revision. Raphaël Galicher (from Soummer, in IDL) -Revision 2.0 2012-04-12 P. Baudoz (IDL version): added pup offset -Revision 3.0 2020-03-10 J. Mazoyer (to python). Replace the MFT with no input offset option -Revision 4.0 2020-04-20 J. Mazoyer. change the normalization. Change dim_pup name to be more

coherent. Made better parameter format check

-Revision 5.0 2022-03-09 J. Mazoyer. 1/2 pixel error in xx0, xx1, uu0 and uu1. Now MFT of clear

pup if fully real.

-Revision 6.0 2022-10-11 J. Mazoyer. Introduced the option to do only the measurement of AA and BB

and the option to do only the matrix multiplication. Matrix multiplication itself is be done separately in mat_mult_mft which allowed GPU maybe. I tried using numba here to save so time but no improvement

Parameters:
image2D array

Entrance image (entrance size in x and y can be different)

real_dim_inputint or tuple of ints of dim 2, default 4

Diameter of the support in pup (can differ from image.shape) Example: real_dim_input = diameter of the pupil in pixel for a padded pupil

dim_outputint or tuple of int of dim 2, default 4

Dimension of the output in pixels (square if int, rectangular if (int, int).

nbres: float or tuple of float of dim 2, default 1

Number of spatial resolution elements (same in both directions if float).

inversebool, default False

Direction of the MFT. If inverse=False, direct mft (default value). If inverse=True, indirect mft.

normstring default ‘backward’

‘backward’, ‘forward’ or ‘ortho’. this is the same paramter as in numpy.fft functions https://numpy.org/doc/stable/reference/routines.fft.html#module-numpy.fft if ‘backward’ no normalisation is done on MFT(inverse = False) and normalisation 1/N is done in MFT(inverse = True) if ‘forward’ 1/N normalisation is done on MFT(inverse = False) and no normalisation is done in MFT(inverse = True) if ‘ortho’ 1/sqrt(N) normalisation is done in both directions. Note that norm = ‘ortho’ allows you to conserve energy between a focal plane and pupil plane The default is ‘backward’ to be consistent with numpy.fft.fft2 and numpy.fft.ifft2

X_offset_inputfloat default 0

position of the 0-frequency pixel in x for the entrance image with respect to the center of the entrance image (real position of the 0-frequency pixel on dim_input_x/2+x0)

Y_offset_inputfloat, default 0

position of the 0-frequency pixel in Y for the entrance image with respect to the center of the entrance image (real position of the 0-frequency pixel on dim_input_y/2+y0)

X_offset_outputfloat, default 0

position of the 0-frequency pixel in x for the output image with respect to the center of the output image (real position of the 0-frequency pixel on dim_output_x/2+x1)

Y_offset_outputfloat, default 0

position of the 0-frequency pixel in Y for the output image with respect to the center of the output image (real position of the 0-frequency pixel on dim_output_y/2+y1)

only_mat_mult: boolean,, default False
if True, we only do the matrix multiplication, but we need AA, BB and norm0 to be provided.

in that case all other parameters are not used. Careful, in this mode, it is assumed that the user is ‘expert’ and no specific error message will be thrown if parameters are wrong. e.g. it will crash if image, AA and BB dimensions are not compatibles

if False : classical MFT, AA, BB and norm0 parameters are not used

AA: complex numpy array, default None

Matrix multiplied in norm0 * ((AA @ image) @ BB). This parameter is only used if only_mat_mult = True

BB: complex numpy array, default None

Matrix multiplied in norm0 * ((AA @ image) @ BB). This parameter is only used if only_mat_mult = True

norm0: float, default None

Normalization value in matrix multiplication norm0 * ((AA @ image) @ BB). This parameter is only used if only_mat_mult = True

returnAABB: boolean, default False

if False, the normal MFT(image is returned) if True, we return AA, BB, norm0 used to do norm0 * ((AA @ image) @ BB).

dtype_complexstring, default ‘complex128’

bit number for the complex arrays in the MFT matrices. Can be ‘complex128’ or ‘complex64’. The latter increases the speed of the mft but at the cost of lower precision.

Returns:
if returnAABB is False:

MFT of the image : complex 2D array.

else:
AA, BB, norm0complex 2D array, complex 2D array, float

terms used in MFT matrix multiplication norm0 * ((AA @ image) @ BB).

Asterix.optics.propagation_functions.mat_mult_mft(image, AA, BB, norm0)

Perform only the Matrix multiplication for the MFT.

It is be done separately in to allow this single line function to be sped up. I tried using the numba compiler on this function (https://numba.pydata.org/) to optimize it, but no improvements. This can probably be optimized with GPU here.

AUTHORJohan Mazoyer

2022-10-11 Creation from MFT

Parameters:
image2D numpy array (complex)

Entrance image

AA2D numpy array (complex)

Matrix multiplied in norm0 * ((AA @ image) @ BB).

BB2D numpy array (complex)

Matrix multiplied in norm0 * ((AA @ image) @ BB).

norm0float

Normalization value in matrix multiplication norm0 * ((AA @ image) @ BB).

Returns:
norm0*((AA@image)@BB)2D numpy array (complex)
Asterix.optics.propagation_functions.prop_fresnel(pup, lam, z, rad, prad, retscale=0, dtype_complex='complex128')

Fresnel propagation of electric field along a distance z in a collimated beam and in free space.

AUTHOR : Raphael Galicher

REVISION HISTORY : - Revision 1.1 2020-01-22 Raphael Galicher Initial revision

Parameters:
pup2D array (complex or real)
if retscale == 0

electric field at z=0 CAUTION : pup has to be centered on (dimpup/2+1,dimpup/2+1) where ‘dimpup’ is the pup array dimension.

else:

Dim of the input array that will be used for pup.

lamfloat

Wavelength in meter.

zfloat

distance of propagation

radfloat

if z>0: entrance beam radius in meter if z<0: output beam radius in meter

pradfloat

if z>0: entrance beam radius in pixel if z<0: output beam radius in pixel

retscaleint 0 or 1:

if not 0, the function returns the scales of the input and output arrays if 0, the function returns the output electric field (see Returns)

dtype_complexstring, default ‘complex128’

bit number for the complex arrays in the exponential. Can be ‘complex128’ or ‘complex64’. The latter increases the speed of the exp but at the cost of lower precision.

Returns:
if retscale is 0:
pup_z2D array (complex)

electric field after propagating in free space along a distance z

dxoutfloat

lateral sampling in the output array

else:
dxfloat

lateral sampling in the input array

dxoutfloat

lateral sampling in the output array

Asterix.optics.propagation_functions.prop_angular_spectrum(pup, lam, z, rad, prad, gamma=2, dtype_complex='complex128')

Angular spectrum propagation of electric field along a distance z in a collimated beam and in free space in close field (small z).

AUTHOR : Johan Mazoyer

REVISION HISTORY : - Revision 1.0 2022-02-15 Johan Mazoyer Initial revision

Parameters:
pup2D array (complex or real)

electric field at z=0 CAUTION : pup has to be centered on (dimpup/2+1,dimpup/2+1) where ‘dimpup’ is the pup array dimension

lamfloat

wavelength in meter

zfloat

distance of propagation in meter

radfloat

entrance beam radius in meter

pradfloat

entrance beam radius in pixel

gammaint >=2

factor of oversizing in the fourier plane in diameter of the pupil (gamma*2*prad is the output dim) optional: default = 2

dtype_complexstring, default ‘complex128’

bit number for the complex arrays in the exponential. Can be ‘complex128’ or ‘complex64’. The latter increases the speed of the exp but at the cost of lower precision.

Returns:
pup_z2D array (complex) of size [2*gamma*prad,2*gamma*prad]

electric field after propagating in free space along a distance z

Asterix.optics.propagation_functions.fft_choosecenter(image, inverse=False, center_pos='bb', norm='backward', dtype_complex='complex128')

FFT Computation. IDL “FFT” routine uses coordinates origin at pixel.

This function uses a coordinate origin at any pixel [k,l], thanks to multiplication by adequate array before using numpy routine “FFT”. Keywords allow convenient origins either at central pixel (‘p’) or between the 4 central pixels (‘b’)

AUTHORS: L.Mugnier, M.Kourdourli, J. Mazoyer

07/09/2022Introduction in asterix (Kourdourli’s version. Based on fftshift2.pro from ONERA’s

IDL library by Laurent Mugnier

07/09/2022 : works for non square array / non even dimensions array Mazoyer

Parameters:
input2D numpy array

Initial array.

inversebool (optional, default False)

Direction of the FFT, inverse == False for direct FFT, inverse == True for inverse FFT.

center_posstring (optional, default ‘bb’)

Option for the origin. Shorthand for specifying the origin center in direct and fourier spaces when manipulating centered arrays.

Direct space Fourier space

pp Central pix Central pix pb Central pix Between 4 central pix bp Between 4 central pix Central pix bb Between 4 central pix Between 4 central pix if dim_i (i = x or y) is even or odd :

Central pix = dim_i // 2 Between 4 central pix: between dim_i // 2 - 1 and dim_i // 2

with // the euclidian division.

normstring default ‘backward’

‘backward’, ‘forward’ or ‘ortho’. this is the same paramter as in numpy.fft functions https://numpy.org/doc/stable/reference/routines.fft.html#module-numpy.fft if ‘backward’ no normalisation is done on MFT(inverse = False) and normalisation 1/N is done in MFT(inverse = True) if ‘forward’ 1/N normalisation is done on MFT(inverse = False) and no normalisation is done in MFT(inverse = True) if ‘ortho’ 1/sqrt(N) normalisation is done in both directions. Note that norm = ‘ortho’ allows you to conserve energy between a focal plane and pupil plane The default is ‘backward’ to be consistent with numpy.fft.fft2 and numpy.fft.ifft2

dtype_complexstring, default ‘complex128’

bit number for the complex arrays in the exponential. Can be ‘complex128’ or ‘complex64’. The latter increases the speed of the exp but at the cost of lower precision. Because numpy fft does not have a dtype parameter, the difference in speed is probably minimum for this function.

Returns:
FFT_array2D numpy array

FFT of input array with respect to the input centering parameters.

Phase and amplitude functions

Asterix.optics.phase_amplitude_functions.roundpupil(dim_pp, prad, grey_pup_bin_factor=1, center_pos='b')

Create a circular pupil.

With grey_pup_bin_factor >1, this creates an oversized pupil that is then rescaled using rebin to dim_pp

AUTHORS : Axel Pottier, Johan Mazoyer 7/9/22 Modified by J Mazoyer to remove the pixel crenellation with rebin and add a better center option 2/10/22 Modified by J Mazoyer to enhance pixel grey_pup_bin_factor factor

Parameters:
dim_ppint

Size of the image (in pixels)

pradfloat

Pupil radius within the image array (in pixels)

grey_pup_bin_factorint (default, 1)

If grey_pup_bin_factor > 1, the pupil is first defined at a very large scale (prad = grey_pup_bin_factor*prad) and then rebinned to the given parameter ‘prad’. This limits the pixel crenellation in the pupil for small pupils. If this option is activated (grey_pup_bin_factor>1) the pupil has to be perfectly centered on the array because binning while keeping the centering is tricky:

-if center_pos is ‘p’, dimpp and grey_pup_bin_factor must both be odd -if center_pos is ‘b’, dimpp and grey_pup_bin_factor must both be even

center_posstring (optional, default ‘b’)

Option for the center pixel. If ‘p’, center on the pixel dim_pp//2. If ‘b’, center in between pixels dim_pp//2 -1 and dim_pp//2, for ‘dim_pp’ odd or even.

Returns:
pupilnormal2D array

Output circular pupil

Asterix.optics.phase_amplitude_functions.shift_phase_ramp(dim_pp, shift_x, shift_y)

Create a phase ramp of size (dim_pp,dim_pp) that can be used as follow to shift one image by (a,b) pixels : shift_im = real(fft(ifft(im)*exp(i phase ramp)))

AUTHOR: Axel Potier

Parameters:
dim_ppint

Size of the phase ramp (in pixels)

shift_xfloat

Shift desired in the x direction (in pixels)

shift_yfloat

Shift desired in the y direction (in pixels)

Returns:
masktot2D array

Phase ramp

Asterix.optics.phase_amplitude_functions.random_phase_map(pupil_rad, dim_image, phaserms, rhoc, slope)

Create a random phase map, whose PSD decreases as f^(-slope).

The average is null and the standard deviation is ‘phaserms’.

AUTHOR: Axel Potier

Parameters:
pupil_rad: float

Radius of the pupil on which the phase rms will be measured.

dim_image: int

Size of the output (can be different from 2*pupil_rad).

phasermsfloat

Standard deviation of the aberration.

rhocfloat

See Borde et Traub 2006.

slopefloat

Slope of the PSD. See Borde et Traub 2006.

Returns:
phase2D array

Static random phase map (or OPD)

Asterix.optics.phase_amplitude_functions.sine_cosine_basis(Nact1D)

For a given number of actuator across the DM, create coefficients for the sin/cos basis.

AUTHOR: Johan Mazoyer

Parameters:
Nact1Dfloat

Number of actuators of a square DM in one of the principal directions.

Returns:
SinCos3D array

Coefficient to apply to DMs to obtain sine and cosine phases. size :[(Nact1D)^2,Nact1D,Nact1D] if even size :[(Nact1D)^2 -1 ,Nact1D,Nact1D] if odd (to remove piston)

Asterix.optics.phase_amplitude_functions.make_apodizer(dim_pp, prad, apodizer_profile, grey_pup_bin_factor=1, center_pos='b')

Return a generic apodizer, by apllying a given transmission profil on the round pupil. The transmission profile must be a fonction of the radial coordinate. The radial coordinate is 0 at the center and 1 on the outer edge, namely at prad.

AUTHOR : Charles Goulas

Parameters:
dim_ppint

Size of the image (in pixels)

pradfloat

Pupil radius within the image array (in pixels)

apodizer_profilefunction.

Apodizer radial transmission. It is a function of the radial coordinate.

grey_pup_bin_factorint (default, 1)

If grey_pup_bin_factor > 1, the pupil is first defined at a very large scale (prad = grey_pup_bin_factor*prad) and then rebinned to the given parameter ‘prad’. This limits the pixel crenellation in the pupil for small pupils. If this option is activated (grey_pup_bin_factor>1) the pupil has to be perfectly centered on the array because binning while keeping the centering is tricky:

-if center_pos is ‘p’, dimpp and grey_pup_bin_factor must both be odd -if center_pos is ‘b’, dimpp and grey_pup_bin_factor must both be even

center_posstring (optional, default ‘b’)

Option for the center pixel. If ‘p’, center on the pixel dim_pp//2. If ‘b’, center in between pixels dim_pp//2 -1 and dim_pp//2, for ‘dim_pp’ odd or even.

Returns:
apodizer_pupil2D array.

Apodizer pupil

Asterix.optics.phase_amplitude_functions.make_spider(dim_pp, starting_point, finishing_point, w_spiders, center_pos='b')

Make a unique spider from starting_point to finishing point, of width w_spiders

AUTHORS : Johan Mazoyer, heavily inspired by Emiel Por in HCIpy [Por2018].

[Por2018]

Por et al. 2018, “High Contrast Imaging for Python (HCIPy):

an open-source adaptive optics and coronagraph simulator”

Parameters:
dim_ppint

Size of the pupil plane (in pixels)

starting_point: tupple of float

(x0,y0) the starting point of the spider in pixel

finishing_point: tupple of float

(x0,y0) the finishing point of the spider in pixel

w_spiders: float

width of the spider in pixel

center_posstring (optional, default ‘b’)

Option for the center pixel. If ‘p’, center on the pixel dim_pp//2. If ‘b’, center in between pixels dim_pp//2 -1 and dim_pp//2, for ‘dim_pp’ odd or even.

Returns:
spider_map2D bool array

spider boolean array

Asterix.optics.phase_amplitude_functions.make_VLT_pup(dim_pp, prad, pupangle=0, spiders=True, grey_pup_bin_factor=1, center_pos='b', reduce_outer_radius=0, add_central_obs=0, add_spider_thickness=0)

Return VLT pup, heavily inspired by HCIpy.

AUTHORS : Johan Mazoyer, heavily inspired by Emiel Por in HCIpy and with help from C. Goulas I used the number in a slide given by Anthony available here: https://www.dropbox.com/s/so0wpq58wh5i5o2/pupil_VLT.pdf?dl=1

Parameters:
dim_ppint

Size of the image (in pixels)

pradfloat

Pupil radius within the image array (in pixels)

pupanglefloat

pupil rotation angle in deg

spidersbool, default True

if False, return the VLT pupil without spiders

grey_pup_bin_factorint, default 1

If grey_pup_bin_factor > 1, the pupil is first defined at a very large scale (prad = grey_pup_bin_factor*prad) and then rebinned to the given parameter ‘prad’. This limits the pixel crenellation in the pupil for small pupils. If this option is activated (grey_pup_bin_factor>1) the pupil has to be perfectly centered on the array because binning while keeping the centering is tricky:

-if center_pos is ‘p’, dimpp and grey_pup_bin_factor must both be odd -if center_pos is ‘b’, dimpp and grey_pup_bin_factor must both be even

center_posstring, default ‘b’

Option for the center pixel. If ‘p’, center on the pixel dim_pp//2. If ‘b’, center in between pixels dim_pp//2 -1 and dim_pp//2, for ‘dim_pp’ odd or even.

reduce_outer_radiusfloat, default 0

reduced diameter of outer radius in fraction of the diameter

add_central_obsfloat, default 0

additional diameter of central obstruction in fraction of the diameter

add_spider_thicknessfloat, default 0

additional spiders width in fraction of the diameter

Returns:
VLTpupil2D numpy array

VLT transmission pupil of shape (pupdiam, pupdiam), filled with 0 and 1

Asterix.optics.phase_amplitude_functions.sphere_apodizer_radial_profile(x)

Compute the transmission radial profile of the SPHERE APO1 apodizer. x is the radial coordinate inside the pupil x = 0 at the center of the pupil and x = 1 on the outer edge This profile has been estimated with a five order polynomial fit. Don’t go inside the central obstruction, namely x < 0.14, as the fit is no longer reliable.

AUTHOR : Charles Goulas

Parameters:
xfloat or array

distance to the pupil center, in fraction of the pupil radius

Returns:
profilefloat or array

corresponding transmission

Asterix.optics.phase_amplitude_functions.make_sphere_apodizer(dim_pp, prad, grey_pup_bin_factor=1, center_pos='b')

Return the SPHERE APO1 apodizer pupil.

AUTHOR : Charles Goulas

Parameters:
dim_ppint

Size of the image (in pixels)

pradfloat

Pupil radius within the image array (in pixels)

grey_pup_bin_factorint, default 1

If grey_pup_bin_factor > 1, the pupil is first defined at a very large scale (prad = grey_pup_bin_factor*prad) and then rebinned to the given parameter ‘prad’. This limits the pixel crenellation in the pupil for small pupils. If this option is activated (grey_pup_bin_factor>1) the pupil has to be perfectly centered on the array because binning while keeping the centering is tricky:

-if center_pos is ‘p’, dimpp and grey_pup_bin_factor must both be odd -if center_pos is ‘b’, dimpp and grey_pup_bin_factor must both be even

center_posstring, default ‘b’

Option for the center pixel. If ‘p’, center on the pixel dim_pp//2. If ‘b’, center in between pixels dim_pp//2 -1 and dim_pp//2, for ‘dim_pp’ odd or even.

Returns:
sphere_apodizer2D array

sphere APO1 apodizer pupil

Asterix.optics.phase_amplitude_functions.make_sphere_lyot(dim_pp, prad, pupangle=0, spiders=True, grey_pup_bin_factor=1, center_pos='b')

Return SPHERE Lyot stop aperture

values of additional central obstruction, spiders size and outer edge obstruction have been estimated by eye on the real lyot stop warning : this lyot stop does not feature the dead actuators patches

AUTHORS : Johan Mazoyer from C. Goulas

Parameters:
dim_ppint

Size of the image (in pixels).

pradfloat

Pupil radius within the image array (in pixels). Careful this is not the radius of the lyot, but the pupil associated to this Lyot.

pupanglefloat

Pupil rotation angle in deg.

spidersbool, default True

if False, return the VLT pupil without spiders

grey_pup_bin_factorint, default 1

If grey_pup_bin_factor > 1, the pupil is first defined at a very large scale (prad = grey_pup_bin_factor*prad) and then rebinned to the given parameter ‘prad’. This limits the pixel crenellation in the pupil for small pupils. If this option is activated (grey_pup_bin_factor>1) the pupil has to be perfectly centered on the array because binning while keeping the centering is tricky:

-if center_pos is ‘p’, dimpp and grey_pup_bin_factor must both be odd -if center_pos is ‘b’, dimpp and grey_pup_bin_factor must both be even

center_posstring, default ‘b’

Option for the center pixel. If ‘p’, center on the pixel dim_pp//2. If ‘b’, center in between pixels dim_pp//2 -1 and dim_pp//2, for ‘dim_pp’ odd or even.

Returns:
SPHERELyotStop2D numpy array

SPHERE Lyot Stop transmission pupil of shape (dim_pp, dim_pp), filled with 0 and 1

Asterix.optics.phase_amplitude_functions.butterworth_circle(dim, size_filter, order=5, xshift=0, yshift=0)

Return a circular Butterworth filter.

AUTHOR: Raphaël Galicher (in IDL)

ILa (to Python)

Parameters:
dimint

Dimension of 2D output array in pixels. If even, filter will be centered on a pixel, but can be shifted to between pixels by using xshift=-0.5 and yshift=-0.5. If uneven, filter will be centered between pixels.

size_filterint

Inverse size of the filter.

orderint

Order of the filter.

xshiftfloat

Shift of filter with respect to its array in the x direction, in pixels.

yshiftfloat

Shift of filter with respect to its array in the y direction, in pixels.

Returns:
butterworth2D array