COMPASS  5.4.4
End-to-end AO simulation tool using GPU acceleration
stellarCoronagraph.py
1 
37 import numpy as np
38 import shesha.config as conf
39 import shesha.constants as scons
41 from shesha.init.coronagraph_init import init_coronagraph, init_mft, mft_multiplication
42 from shesha.supervisor.components.targetCompass import TargetCompass
43 from sutraWrap import StellarCoronagraph
44 from carmaWrap import context
45 class StellarCoronagraphCompass(GenericCoronagraph):
46  """ Class supervising stellar coronagraph component
47 
48  Attributes:
49  _spupil: (np.ndarray[ndim=2, dtype=np.float32]): Telescope pupil mask
50 
51  _pupdiam : (int): Number of pixels along the pupil diameter
52 
53  _dim_image :(int): Coronagraphic image dimension
54 
55  _p_corono: (Param_corono): Coronagraph parameters
56 
57  _target: (TargetCompass): Compass Target used as input for the coronagraph
58 
59  _norm_img : (float): Normalization factor for coronagraphic image
60 
61  _norm_psf : (float): Normalization factor for PSF
62 
63  _coronagraph: (SutraCoronagraph): Sutra coronagraph instance
64 
65  _wav_vec: (np.ndarray[ndim=1, dtype=np.float32]): Vector of wavelength
66 
67  _AA_apod_to_fpm: (np.ndarray[ndim=3, dtype=np.complex64]): MFT matrix for focal plane
68 
69  _BB_apod_to_fpm: (np.ndarray[ndim=2, dtype=np.complex64]): MFT matrix for focal plane
70 
71  _norm0_apod_to_fpm: (np.ndarray[ndim=2, dtype=np.complex64]): MFT matrix for focal plane
72 
73  _AA_fpm_to_lyot: (np.ndarray[ndim=3, dtype=np.complex64]): MFT matrix for lyot plane
74 
75  _BB_fpm_to_lyot: (np.ndarray[ndim=2, dtype=np.complex64]): MFT matrix for lyot plane
76 
77  _norm0_fpm_to_lyot: (np.ndarray[ndim=2, dtype=np.complex64]): MFT matrix for lyot plane
78 
79  _AA_lyot_to_image_c: (np.ndarray[ndim=2, dtype=np.complex64]): MFT matrix for image computation (centered on pixel)
80 
81  _BB_lyot_to_image_c: (np.ndarray[ndim=2, dtype=np.complex64]): MFT matrix for psf computation (centered on pixel)
82 
83  _norm0_lyot_to_image_c: (np.ndarray[ndim=2, dtype=np.complex64]): MFT matrix for psf computation (centered on pixel)
84 
85  _AA_lyot_to_image: (np.ndarray[ndim=2, dtype=np.complex64]): MFT matrix for image computation
86 
87  _BB_lyot_to_image: (np.ndarray[ndim=2, dtype=np.complex64]): MFT matrix for psf computation
88 
89  _norm0_lyot_to_image: (np.ndarray[ndim=2, dtype=np.complex64]): MFT matrix for psf computation
90 
91  _indices_pup: (tuple): Tuple of ndarray containing X and Y indices of illuminated
92  pixels in the pupil
93  """
94  def __init__(self, context: context, targetCompass: TargetCompass, p_corono: conf.Param_corono,
95  p_geom: conf.Param_geom):
96  """ Initialize a stellar coronagraph instance
97 
98  Args:
99  context: (CarmaWrap.context): GPU context
100 
101  targetCompass: (TargetCompass): Compass Target used as input for the coronagraph
102 
103  p_corono: (Param_corono): Coronagraph parameters
104 
105  p_geom: (Param_geom): Compass geometry parameters
106  """
107 
108  init_coronagraph(p_corono, p_geom.pupdiam)
109  GenericCoronagraph.__init__(self, p_corono, p_geom, targetCompass)
110  self._wav_vec_wav_vec = p_corono._wav_vec
111 
112  self._AA_apod_to_fpm, self._BB_apod_to_fpm, self._norm0_apod_to_fpm_norm0_apod_to_fpm = init_mft(self._p_corono_p_corono,
113  self._pupdiam_pupdiam,
114  planes='apod_to_fpm')
115  self._AA_fpm_to_lyot, self._BB_fpm_to_lyot, self._norm0_fpm_to_lyot_norm0_fpm_to_lyot = init_mft(self._p_corono_p_corono,
116  self._pupdiam_pupdiam,
117  planes='fpm_to_lyot')
118  self._AA_lyot_to_image, self._BB_lyot_to_image, self._norm0_lyot_to_image_norm0_lyot_to_image = init_mft(self._p_corono_p_corono,
119  self._pupdiam_pupdiam,
120  planes='lyot_to_image')
121  self._AA_lyot_to_image_c, self._BB_lyot_to_image_c, self._norm0_lyot_to_image_c_norm0_lyot_to_image_c = init_mft(self._p_corono_p_corono,
122  self._pupdiam_pupdiam,
123  planes='lyot_to_image',
124  center_on_pixel=True)
125  self._coronagraph_coronagraph_coronagraph = StellarCoronagraph(context, self._target_target.sources[0],
126  self._dim_image_dim_image, self._dim_image_dim_image,
127  self._p_corono_p_corono._dim_fpm, self._p_corono_p_corono._dim_fpm,
128  self._wav_vec_wav_vec, self._wav_vec_wav_vec.size,
129  self._p_corono_p_corono._babinet_trick, 0)
130 
131  self._coronagraph_coronagraph_coronagraph.set_mft(self._AA_lyot_to_image,
132  self._BB_lyot_to_image,
133  self._norm0_lyot_to_image_norm0_lyot_to_image, scons.MftType.IMG)
134  self._coronagraph_coronagraph_coronagraph.set_mft(self._AA_lyot_to_image_c,
135  self._BB_lyot_to_image_c,
136  self._norm0_lyot_to_image_c_norm0_lyot_to_image_c, scons.MftType.PSF)
137  self._coronagraph_coronagraph_coronagraph.set_mft(self._AA_apod_to_fpm,
138  self._BB_apod_to_fpm,
139  self._norm0_apod_to_fpm_norm0_apod_to_fpm, scons.MftType.FPM)
140  self._coronagraph_coronagraph_coronagraph.set_mft(self._AA_fpm_to_lyot,
141  self._BB_fpm_to_lyot,
142  self._norm0_fpm_to_lyot_norm0_fpm_to_lyot, scons.MftType.LYOT)
143 
144  self._coronagraph_coronagraph_coronagraph.set_apodizer(self._p_corono_p_corono._apodizer)
145  self._coronagraph_coronagraph_coronagraph.set_lyot_stop(self._p_corono_p_corono._lyot_stop)
146  fpm = np.rollaxis(np.array(self._p_corono_p_corono._focal_plane_mask), 0, 3)
147  if self._p_corono_p_corono._babinet_trick:
148  fpm = 1. - fpm
149  self._coronagraph_coronagraph_coronagraph.set_focal_plane_mask(fpm)
150  self._compute_normalization_compute_normalization()
151 
152  def _compute_normalization(self):
153  """ Compute the normalization factor of coronagraphic images
154  """
155  self._target_target.reset_tar_phase(0)
156  self._coronagraph_coronagraph_coronagraph.compute_image_normalization()
157  self._norm_img_norm_img_norm_img = np.max(self.get_imageget_image(expo_type=scons.ExposureType.SE))
158  self.compute_psfcompute_psf(accumulate=False)
159  self._norm_psf_norm_psf_norm_psf = np.max(self.get_psfget_psf(expo_type=scons.ExposureType.SE))
def get_psf(self, *str expo_type=scons.ExposureType.LE)
Return the psf.
def compute_psf(self, *bool accumulate=True)
Compute the SE psf, and accumulate it in the LE image.
def get_image(self, *str expo_type=scons.ExposureType.LE)
Return the coronagraphic image.
def __init__(self, context context, TargetCompass targetCompass, conf.Param_corono p_corono, conf.Param_geom p_geom)
Parameter classes for COMPASS.
Numerical constants for shesha and config enumerations for safe-typing.
Definition: constants.py:1
Initialization of a Coronagraph object.