Source code for webbpsf.wfirst

"""
==============================
WFIRST Instruments (pre-alpha)
==============================

WARNING: This model has not yet been validated against other PSF
         simulations, and uses several approximations (e.g. for
         mirror polishing errors, which are taken from HST).
"""

import os.path
import poppy
import numpy as np
from . import webbpsf_core
from scipy.interpolate import griddata
from astropy.io import fits
import astropy.units as u
import logging

_log = logging.getLogger('webbpsf')
import pprint


class WavelengthDependenceInterpolator(object):
    """WavelengthDependenceInterpolator can be configured with
    `n_zernikes` worth of Zernike coefficients at up to `n_wavelengths`
    wavelengths, and will let you `get_aberration_terms` for any
    wavelength in range interpolated linearly between measured/known
    points
    """

    def __init__(self, n_wavelengths=16, n_zernikes=22):
        self._n_wavelengths = n_wavelengths
        self._n_zernikes = n_zernikes
        self._aberration_terms = np.zeros((n_wavelengths, n_zernikes), dtype=np.float64)
        self._wavelengths = []

    def set_aberration_terms(self, wavelength, zernike_array):
        """Supply a reference `wavelength` and a `zernike_array`
        (of length `n_zernikes`) where the aberration is known
        """
        n_wavelengths_set = len(self._wavelengths)
        if wavelength not in self._wavelengths and n_wavelengths_set < self._n_wavelengths:
            self._wavelengths.append(wavelength)
            aberration_row_idx = n_wavelengths_set  # which is now index of last row
        elif wavelength in self._wavelengths:
            aberration_row_idx = self._wavelengths.index(wavelength)
        else:
            # can't add more wavelengths without allocating new _aberration_terms array
            raise ValueError("Already have information at {} wavelengths "
                             "(pass larger n_wavelengths to __init__?)".format(self._n_wavelengths))
        if len(zernike_array) != self._n_zernikes:
            raise ValueError("Expected {} aberration terms (pass different "
                             "n_zernikes to __init__?)".format(self._n_zernikes))
        self._aberration_terms[aberration_row_idx] = zernike_array

    def get_aberration_terms(self, wavelength):
        """Return the Zernike coefficients as interpolated for this
        `wavelength`"""
        # return array of length n_zernikes interpolated for this wavelength
        if wavelength in self._wavelengths:
            # aberration known exactly for this wavelength
            aberration_row_idx = self._wavelengths.index(wavelength)
            return self._aberration_terms[aberration_row_idx]
        else:
            # we have to interpolate @ this wavelength
            aberration_terms = griddata(self._wavelengths, self._aberration_terms, wavelength, method='linear')
            if np.any(np.isnan(aberration_terms)):
                if isinstance(wavelength, u.Quantity):
                    wavelength = wavelength.to(u.m).value
                wavelength_closest = np.clip(wavelength, np.min(self._wavelengths), np.max(self._wavelengths))
                _log.warn("Attempted to get aberrations at wavelength {:.2g} "
                          "outside the range of the reference data; clipping to closest wavelength {:.2g}".format(
                    wavelength, wavelength_closest))

                aberration_terms = griddata(self._wavelengths, self._aberration_terms, wavelength_closest,
                                            method='linear')
            return aberration_terms

class FieldDependentAberration(poppy.ZernikeWFE):
    """FieldDependentAberration incorporates aberrations that
    are interpolated in wavelength, x, and y pixel positions by
    computing the Zernike coefficients for a particular wavelength
    and position.
    """

    """By default, `get_aberration_terms` will zero out Z1, Z2, and Z3
    (piston, tip, and tilt) as they are not meaningful for telescope
    PSF calculations (the former is irrelevant, the latter two would
    be handled by a distortion solution). Change
    `_omit_piston_tip_tilt` to False to include the Z1-3 terms."""
    _omit_piston_tip_tilt = True
    _field_position = None

    def __init__(self, pixel_width, pixel_height,
                 name="Field-dependent Aberration", radius=1.0, oversample=1, interp_order=3):
        self.pixel_width, self.pixel_height = pixel_width, pixel_height
        self.field_position = pixel_width // 2, pixel_height // 2
        self._wavelength_interpolators = {}
        self.pupil_diam = radius * 2.0
        super(FieldDependentAberration, self).__init__(
            name=name,
            verbose=True,
            radius=radius,
            oversample=oversample,
            interp_order=interp_order
        )

    def get_opd(self, wave):
        """Set the Zernike coefficients (for ZernikeWFE.getOPD) based
        on the wavelength of the incoming wavefront and the pixel
        position
        """
        if not isinstance(wave, poppy.Wavefront):
            wavelength = wave
        else:
            wavelength = wave.wavelength
        self.coefficients = wavelength * self.get_aberration_terms(wavelength)
        return super(FieldDependentAberration, self).get_opd(wave)

    @property
    def field_position(self):
        return self._field_position

    @field_position.setter
    def field_position(self, position):
        """Set the x and y pixel position on the detector for which to
        interpolate aberrations"""
        x_pixel, y_pixel = position
        if x_pixel > self.pixel_width or x_pixel < 0:
            raise ValueError("Requested pixel_x position lies outside "
                             "the detector width ({})".format(x_pixel))
        if y_pixel > self.pixel_height or y_pixel < 0:
            raise ValueError("Requested pixel_y position lies outside "
                             "the detector height ({})".format(y_pixel))

        self._field_position = x_pixel, y_pixel

    def add_field_point(self, x_pixel, y_pixel, interpolator):
        """Supply a wavelength-space interpolator for a pixel position
        on the detector"""
        self._wavelength_interpolators[(x_pixel, y_pixel)] = interpolator

    def get_aberration_terms(self, wavelength):
        """Supply the Zernike coefficients for the aberration based on
        the wavelength and pixel position on the detector"""
        if self.field_position in self._wavelength_interpolators:
            # short path: this is a known point
            interpolator = self._wavelength_interpolators[self.field_position]
            coefficients = interpolator.get_aberration_terms(wavelength)
        else:
            # get aberrations at all field points
            field_points, aberration_terms = [], []
            for field_point_coords, point_interpolator in self._wavelength_interpolators.items():
                field_points.append(field_point_coords)
                aberration_terms.append(point_interpolator.get_aberration_terms(wavelength))
            aberration_array = np.asarray(aberration_terms)
            assert len(aberration_array.shape) == 2, "computed aberration array is not 2D " \
                                                     "(inconsistent number of Zernike terms " \
                                                     "at each point?)"

            field_position = tuple(np.clip(self.field_position, 4, 4092))
            coefficients = griddata(
                np.asarray(field_points),
                np.asarray(aberration_terms),
                field_position,
                method='linear'
            )
            if np.any(np.isnan(coefficients)):
                raise RuntimeError("Could not get aberrations for input field point")
        if self._omit_piston_tip_tilt:
            _log.debug("Omitting piston/tip/tilt")
            coefficients[:3] = 0.0  # omit piston, tip, and tilt Zernikes
        return coefficients


def _load_wfi_detector_aberrations(filename):
    from astropy.io import ascii
    zernike_table = ascii.read(filename)
    detectors = {}

    def build_detector_from_table(number, zernike_table):
        """Build a FieldDependentAberration optic for a detector using
        Zernikes Z1-Z22 at various wavelengths and field points"""
        single_detector_info = zernike_table[zernike_table['sca'] == number]
        field_points = set(single_detector_info['field_point'])
        interpolators = {}
        detector = FieldDependentAberration(
            4096,
            4096,
            radius=WFIRSTInstrument.PUPIL_RADIUS,
            name="Field Dependent Aberration (SCA{:02})".format(number)
        )
        for field_id in field_points:
            field_point_rows = single_detector_info[single_detector_info['field_point'] == field_id]
            local_x, local_y = field_point_rows[0]['local_x'], field_point_rows[0]['local_y']
            interpolator = build_wavelength_dependence(field_point_rows)

            midpoint_pixel = 4096 / 2
            # (local_x in mm / 10 um pixel size) -> * 1e2
            # local_x and _y range from -20.44 to +20.44, so adding to the midpoint pixel
            # makes sense to place (-20.44, -20.44) at (4, 4)
            pixx, pixy = (round(midpoint_pixel + local_x * 1e2),
                          round(midpoint_pixel + local_y * 1e2))

            detector.add_field_point(pixx, pixy, interpolator)
        return detector

    def build_wavelength_dependence(rows):
        """Build an interpolator object that interpolates Z1-Z22 in
        wavelength space"""
        wavelengths = set(rows['wavelength'])
        interpolator = WavelengthDependenceInterpolator(n_wavelengths=len(wavelengths),
                                                        n_zernikes=22)
        for row in rows:
            z = np.zeros(22)
            for idx in range(22):
                z[idx] = row['Z{}'.format(idx + 1)]
            interpolator.set_aberration_terms(row['wavelength'] * 1e-6, z)

        return interpolator

    detector_ids = set(zernike_table['sca'])
    for detid in detector_ids:
        detectors["SCA{:02}".format(detid)] = build_detector_from_table(detid, zernike_table)

    return detectors


class WFIRSTInstrument(webbpsf_core.SpaceTelescopeInstrument):
    PUPIL_RADIUS = 2.4 / 2.0
    """
    WFIRSTInstrument contains data and functionality common to WFIRST
    instruments, such as setting the pupil shape
    """
    telescope = "WFIRST"

    def __init__(self, *args, **kwargs):
        super(WFIRSTInstrument, self).__init__(*args, **kwargs)

    # slightly different versions of the following two functions
    # from the parent superclass
    # in order to interface with the FieldDependentAberration class
    @property
    def detector_position(self):
        """The pixel position in (X, Y) on the detector"""
        return self._detectors[self._detector].field_position

    @detector_position.setter
    def detector_position(self, position):
        # exact copy of superclass function except we save the
        # into a different location.
        try:
            x, y = map(int, position)
        except ValueError:
            raise ValueError("Detector pixel coordinates must be pairs of nonnegative numbers, "
                             "not {}".format(position))
        if x < 0 or y < 0:
            raise ValueError("Detector pixel coordinates must be nonnegative integers")
        if x > self._detector_npixels - 1 or y > self._detector_npixels - 1:
            raise ValueError("The maximum allowed detector pixel "
                             "coordinate value is {}".format(self._detector_npixels - 1))

        self._detectors[self._detector].field_position = (int(position[0]), int(position[1]))

    def _get_aberrations(self):
        """Get the OpticalElement that applies the field-dependent
        optical aberrations. (Called in _get_optical_system.)"""
        return self._detectors[self._detector]

    def _get_fits_header(self, result, options):
        """Populate FITS Header keywords"""
        super(WFIRSTInstrument, self)._get_fits_header(result, options)
        result[0].header['DETXPIXL'] = (self.detector_position[0],
                                        'X pixel position (for field dependent aberrations)')
        result[0].header['DETYPIXL'] = (self.detector_position[1],
                                        'Y pixel position (for field dependent aberrations)')
        result[0].header['DETECTOR'] = (self.detector, 'Detector selected')


class WFIPupilController:
    """
    This is a helper class for the WFI and is used to swap in
    the correct pupil each time the detector is changed.
    The pupil depends on the selected detector, filter and the
    pupil_mask flag provided by the user.
    The user should not interact with this class directly, only
    through the API provided through the WFI class.
    """

    def __init__(self):
        self._datapath = None
        self._pupil_basepath = None

        self._pupil = None

        # Paths to the two possible pupils. The correct one is selected based on requested
        # wavelengths in _validate_config()
        self._unmasked_pupil_path = None
        self._masked_pupil_path = None

        # List of filters that need the masked pupil
        self._masked_filters = ['F184']

        # Flag to en-/disable automatic selection of the appropriate pupil_mask
        self.auto_pupil = True

        self._pupil_mask = "AUTO"
        self.pupil_mask_list = ['AUTO', 'COLD_PUPIL', 'UNMASKED']
        self._currently_masked = False

    def set_base_path(self, datapath):
        """
        Sets the path to the WebbPSF data files.
        This should be set before this class is used.
        Parameters
        ----------
        datapath : string
            Path to WebbPSF-WFI data files
        """
        self._datapath = datapath
        self._pupil_basepath = os.path.join(self._datapath, "pupils")

    @property
    def pupil(self):
        return self._pupil

    @pupil.setter
    def pupil(self, value):
        self._pupil = value

    @property
    def pupil_mask(self):
        """
        pupil_mask types:
        - "AUTO":
            Automatically select pupil
        - "COLD_PUPIL":
            Masked pupil override
        - "UNMASKED":
            Unmasked pupil override
        """
        return self._pupil_mask

    @pupil_mask.setter
    def pupil_mask(self, name):
        """
        Set the pupil mask

        Parameters
        ------------
        name : string
            Name of setting.
            Settings:
                - "AUTO":
                    Automatically select pupil
                - "COLD_PUPIL":
                    Masked pupil override
                - "UNMASKED":
                    Unmasked pupil override
        """
        if name and isinstance(name, str):
            name = name.upper()
            if "AUTO" == name:
                self.auto_pupil = True
                _log.info("Using default pupil mask.")
            elif "COLD_PUPIL" == name:
                self.auto_pupil = False
                _log.info("Using custom pupil mask: Masked Pupil.")
            elif "UNMASKED" == name:
                self.auto_pupil = False
                _log.info("Using custom pupil mask: Unmasked Pupil.")
            else:
                raise ValueError("Instrument {0} doesn't have a pupil mask called '{1}'.".format(self.name, name))
        else:
            raise ValueError("Pupil mask setting is not valid or empty.")
        self._pupil_mask = name

        self._update_pupil()

    def update_pupil_path(self, detector):
        """
        Update the masked and unmasked pupil paths according to the SCA selected
        """
        if self._pupil_basepath is None:
            raise Exception("update_pupil_path called before pupil file path is set")
        if detector is None:
            raise ValueError("Detector was not set when trying to set pupil file path")
        if 'SCA' not in detector:
            raise ValueError("Unidentified detector selected, could not assign pupil")

        detector = detector[:3] + str(int((detector[3:])))  # example "SCA01" -> "SCA1"

        self._unmasked_pupil_path = os.path.join(self._pupil_basepath,
                                                 '{}_rim_mask.fits.gz'.format(detector))

        self._masked_pupil_path = os.path.join(self._pupil_basepath,
                                               '{}_full_mask.fits.gz'.format(detector))
        self._update_pupil()

    def _update_pupil(self):
        """
        Update the actual pupil by setting the pupil variable
        to the correct pupil path.
        """
        if self._pupil_basepath is None:
            raise Exception("update pupil called before pupil file path is set")

        if 'AUTO' == self.pupil_mask:
            if self._currently_masked:
                self.pupil = self._masked_pupil_path
            else:
                self.pupil = self._unmasked_pupil_path
        elif "COLD_PUPIL" == self.pupil_mask:
            self.pupil = self._masked_pupil_path
        elif "UNMASKED" == self.pupil_mask:
            self.pupil = self._unmasked_pupil_path
        else:
            raise ValueError("Pupil mask setting is not valid or empty.")

    def validate_pupil(self, filter, **kwargs):
        """Validates that the WFI is configured sensibly

        This mainly consists of selecting the masked or unmasked pupil
        appropriately based on the wavelengths requested.
        """
        if self.auto_pupil:
            if filter in self._masked_filters:
                # use masked pupil optic
                self.pupil = self._masked_pupil_path
                _log.info("Using the masked WFI pupil shape based on filter requested")
            else:
                # use unmasked pupil optic
                self.pupil = self._unmasked_pupil_path
                _log.info("Using the unmasked WFI pupil shape based on filter requested")
        else:
            # If the user has set the pupil to a custom value, let them worry about the
            # correct shape it should have
            pass

    def remove_pupil_mask_override(self):
        _log.info("Removing custom pupil mask")
        self.pupil_mask = 'AUTO'


[docs]class WFI(WFIRSTInstrument): """ WFI represents the WFIRST wide field imager for the WFIRST mission WARNING: This model has not yet been validated against other PSF simulations, and uses several approximations (e.g. for mirror polishing errors, which are taken from HST). """ def __init__(self): """ Initiate WFI Parameters ----------- set_pupil_mask_on : bool or None Set to True or False to force using or not using the cold pupil mask, or to None for the automatic behavior. """ pixelscale = 110e-3 # arcsec/px, WFIRST-AFTA SDT report final version (p. 91) # Initialize the pupil controller self._pupil_controller = WFIPupilController() super(WFI, self).__init__("WFI", pixelscale=pixelscale) self._pupil_controller.set_base_path(self._datapath) self.pupil_mask_list = self._pupil_controller.pupil_mask_list self._detector_npixels = 4096 self._detectors = _load_wfi_detector_aberrations(os.path.join(self._datapath, 'wim_zernikes_cycle8.csv')) assert len(self._detectors.keys()) > 0 self.detector = 'SCA01' self.opd_list = [ os.path.join(self._WebbPSF_basepath, 'upscaled_HST_OPD.fits'), ] self.pupilopd = self.opd_list[-1] def _validate_config(self, **kwargs): """Validates that the WFI is configured sensibly This mainly consists of selecting the masked or unmasked pupil appropriately based on the wavelengths requested. """ self._pupil_controller.validate_pupil(self.filter, **kwargs) super(WFI, self)._validate_config(**kwargs) @WFIRSTInstrument.detector.setter def detector(self, value): if value.upper() not in self.detector_list: raise ValueError("Invalid detector. Valid detector names are: {}".format(', '.join(self.detector_list))) self._detector = value.upper() self._pupil_controller.update_pupil_path(self.detector) @property def pupil(self): return self._pupil_controller.pupil @pupil.setter def pupil(self, value): # self._pupil_controller is not available at initiation thus # we must ignore any assignments at super(WFI, self).__init__(...) if self._pupil_controller: self._pupil_controller.pupil = value @property def pupil_mask(self): return self._pupil_controller.pupil_mask @WFIRSTInstrument.filter.setter def filter(self, value): value = value.upper() # force to uppercase if value not in self.filter_list: raise ValueError("Instrument %s doesn't have a filter called %s." % (self.name, value)) self._filter = value self._pupil_controller.validate_pupil(self.filter) @pupil_mask.setter def pupil_mask(self, name): """ Set the pupil mask Parameters ------------ name : string Name of setting. Settings: - "AUTO": Automatically select pupil - "COLD_PUPIL": Masked pupil override - "UNMASKED": Unmasked pupil override """ self._pupil_controller.pupil_mask = name def _addAdditionalOptics(self, optsys, **kwargs): return optsys, False, None @property def _unmasked_pupil_path(self): return self._pupil_controller._unmasked_pupil_path @property def _masked_pupil_path(self): return self._pupil_controller._masked_pupil_path
[docs]class CGI(WFIRSTInstrument): """ WFIRST Coronagraph Instrument Simulates the PSF of the WFIRST coronagraph. Current functionality is limited to the Shaped Pupil Coronagraph (SPC) observing modes, and these modes are only simulated with static, unaberrated wavefronts, without relay optics and without DM control. The design respresented here is an approximation to a baseline concept, and will be subject to change based on trades studies and technology development. Parameters ---------- mode : str CGI observing mode. If not specified, the __init__ function will set this to a default mode 'CHARSPC_F660' pixelscale : float Detector pixelscale. If not specified, the pixelscale will default to 0.02 arcsec for configurations usint the IMAGER camera and 0.025 arcsec for the IFS. fov_arcsec : float Field of view in arcseconds. If not specified, the field of view will default to 3.20 arcsec for the IMAGER camera and 1.76 arcsec for the IFS. """ camera_list = ['IMAGER', 'IFS'] filter_list = ['F660', 'F721', 'F770', 'F890'] apodizer_list = ['CHARSPC', 'DISKSPC'] fpm_list = ['CHARSPC_F660_BOWTIE', 'CHARSPC_F770_BOWTIE', 'CHARSPC_F890_BOWTIE', 'DISKSPC_F721_ANNULUS'] lyotstop_list = ['LS30D88'] _mode_table = { # MODE CAMERA FILTER APODIZER FPM LYOT STOP 'CHARSPC_F660': ('IFS', 'F660', 'CHARSPC', 'CHARSPC_F660_BOWTIE', 'LS30D88'), 'CHARSPC_F770': ('IFS', 'F770', 'CHARSPC', 'CHARSPC_F770_BOWTIE', 'LS30D88'), 'CHARSPC_F890': ('IFS', 'F890', 'CHARSPC', 'CHARSPC_F890_BOWTIE', 'LS30D88'), 'DISKSPC_F721': ('IMAGER', 'F721', 'DISKSPC', 'DISKSPC_F721_ANNULUS', 'LS30D88')} def __init__(self, mode=None, pixelscale=None, fov_arcsec=None, apply_static_opd=False): super(CGI, self).__init__("CGI", pixelscale=pixelscale) self._detector_npixels = 1024 self._detectors = {camera: 'placeholder' for camera in self.camera_list} self.pupil_mask_list = self.lyotstop_list # alias for use in webbpsf_core self.image_mask_list = self.fpm_list # alias for use in webbpsf_core self.pupil = os.path.join(self._WebbPSF_basepath, 'AFTA_CGI_C5_Pupil_onax_256px_flip.fits') if apply_static_opd: self.pupilopd = os.path.join(self._WebbPSF_basepath, 'CGI', 'OPD', 'CGI_static_OPD.fits') else: self.pupilopd = None self.aberration_optic = None self.options = {'force_coron': True} # Allow the user to pre-emptively override the default instrument FoV and pixel scale if fov_arcsec is not None: self.fov_arcsec = fov_arcsec self._override_fov = True else: self._override_fov = False if pixelscale is not None: self._pixelscale = pixelscale self._override_pixelscale = True else: self._override_pixelscale = False if mode is None: self.print_mode_table() _log.info("Since the mode was not specified at instantiation, defaulting to CHARSPC_F660") self.mode = 'CHARSPC_F660' else: self.mode = mode @property def camera(self): """Currently selected camera name""" return self._camera @camera.setter def camera(self, value): value = value.upper() # force to uppercase if value not in self.camera_list: raise ValueError("Instrument {0} doesn't have a camera called {1}.".format(self.name, value)) self._camera = value if value == 'IMAGER': if not hasattr(self, 'fov_arcsec') or not self._override_fov: self.fov_arcsec = 3.2 if not hasattr(self, 'pixelscale') or not self._override_pixelscale: self.pixelscale = 0.020 # Nyquist at 465 nm else: # default to 'IFS' if not hasattr(self, 'fov_arcsec') or not self._override_fov: self.fov_arcsec = 2 * 0.82 # 2015 SDT report, Section 3.4.1.1.1: # IFS has 76 lenslets across the (2 x 0.82) arcsec FoV. if not hasattr(self, 'pixelscale') or not self._override_pixelscale: self.pixelscale = 0.025 # Nyquist at 600 nm # for CGI, there is one detector per camera and it should be set automatically. @property def detector(self): return self.camera @detector.setter def detector(self, value): raise RuntimeError("Can't set detector directly for CGI; set camera instead.") @property def filter(self): """Currently selected filter name""" return self._filter @filter.setter def filter(self, value): value = value.upper() # force to uppercase if value not in self.filter_list: raise ValueError("Instrument {0} doesn't have a filter called {1}.".format(self.name, value)) self._filter = value @property def apodizer(self): """Currently selected apodizer name""" return self._apodizer @apodizer.setter def apodizer(self, value): value = value.upper() # force to uppercase if value not in self.apodizer_list: raise ValueError("Instrument {0} doesn't have a apodizer called {1}.".format(self.name, value)) self._apodizer = value if value == 'DISKSPC': self._apodizer_fname = \ os.path.join(self._datapath, "optics/DISKSPC_SP_256pix.fits.gz") else: # for now, default to CHARSPC self._apodizer_fname = \ os.path.join(self._datapath, "optics/CHARSPC_SP_256pix.fits.gz") @property def fpm(self): """Currently selected FPM name""" return self._fpm @fpm.setter def fpm(self, value): value = value.upper() # force to uppercase if value not in self.fpm_list: raise ValueError("Instrument {0} doesn't have a FPM called {1}.".format(self.name, value)) self._fpm = value if value.startswith('DISKSPC'): self._fpmres = 3 self._owa = 20. self._Mfpm = int(np.ceil(self._fpmres * self._owa)) self._fpm_fname = \ os.path.join(self._datapath, "optics/DISKSPC_FPM_65WA200_360deg_-_FP1res{0:d}_evensamp_D{1:03d}_{2:s}.fits.gz".format( self._fpmres, 2 * self._Mfpm, self.filter)) else: self._fpmres = 4 self._owa = 9. self._Mfpm = int(np.ceil(self._fpmres * self._owa)) self._fpm_fname = \ os.path.join(self._datapath, "optics/CHARSPC_FPM_25WA90_2x65deg_-_FP1res{0:d}_evensamp_D{1:03d}_{2:s}.fits.gz".format( self._fpmres, 2 * self._Mfpm, self.filter)) @property def lyotstop(self): """Currently selected Lyot stop name""" return self._lyotstop @lyotstop.setter def lyotstop(self, value): # preserve case for this one since we're used to that with the lyot mask names if value not in self.lyotstop_list: raise ValueError("Instrument {0} doesn't have a Lyot mask called {1}.".format(self.name, value)) self._lyotstop = value self._lyotstop_fname = os.path.join(self._datapath, "optics/SPC_LS_30D88_256pix.fits.gz") @property def mode_list(self): """Available Observation Modes""" keys = self._mode_table.keys() keys = sorted(keys) return keys # mode works differently since it's a meta-property that affects the other ones: @property def mode(self): """Currently selected mode name""" for modename, settings in self._mode_table.items(): if (self.camera == settings[0].upper() and self.filter == settings[1].upper() and self.apodizer == settings[2].upper() and self.fpm == settings[3].upper() and self.lyotstop == settings[4]): return modename return 'Custom' @mode.setter def mode(self, value): if value not in self.mode_list: raise ValueError("Instrument {0} doesn't have a mode called {1}.".format(self.name, value)) settings = self._mode_table[value] self.camera = settings[0] self.filter = settings[1] self.apodizer = settings[2] self.fpm = settings[3] self.lyotstop = settings[4] _log.info('Set the following optical configuration:') _log.info('camera = {0}, filter = {1}, apodizer = {2}, fpm = {3}, lyotstop = {4}'.format(\ self.camera, self.filter, self.apodizer, self.fpm, self.lyotstop))
[docs] def print_mode_table(self): """Print the table of observing mode options and their associated optical configuration""" _log.info("Printing the table of WFIRST CGI observing modes supported by WebbPSF.") _log.info("Each is defined by a combo of camera, filter, apodizer, " "focal plane mask (FPM), and Lyot stop settings:") _log.info(pprint.pformat(self._mode_table))
@property def detector_position(self): """The pixel position in (X, Y) on the detector""" return 512, 512 @detector_position.setter def detector_position(self, position): raise RuntimeError("Detector position not adjustable for CGI") def _validate_config(self, **kwargs): super(CGI, self)._validate_config(**kwargs) def _addAdditionalOptics(self, optsys, oversample=4): """Add coronagraphic or spectrographic optics for WFIRST CGI.""" trySAM = False if ('pupil_shift_x' in self.options and self.options['pupil_shift_x'] != 0) or \ ('pupil_shift_y' in self.options and self.options['pupil_shift_y'] != 0): shift = (self.options['pupil_shift_x'], self.options['pupil_shift_y']) else: shift = None # Add the shaped pupil apodizer optsys.add_pupil(transmission=self._apodizer_fname, name=self.apodizer, shift=None) # Add the FPM optsys.add_image(transmission=self._fpm_fname, name=self.fpm) # Add Lyot stop self.pupil_mask = self.lyotstop optsys.add_pupil(transmission=self._lyotstop_fname, name=self.lyotstop, shift=shift) # Cast as MatrixFTCoronagraph; this configures the detector occ_box_size = 1. mft_optsys = poppy.MatrixFTCoronagraph(optsys, oversample=oversample, occulter_box=occ_box_size) return mft_optsys, trySAM, occ_box_size def _get_aberrations(self): """Get the OpticalElement that applies the field-dependent optical aberrations. (Called in _get_optical_system.)""" return None def _get_fits_header(self, result, options): """Populate FITS Header keywords""" super(WFIRSTInstrument, self)._get_fits_header(result, options) pupil_hdr = fits.getheader(self.pupil) apodizer_hdr = fits.getheader(self._apodizer_fname) fpm_hdr = fits.getheader(self._fpm_fname) lyotstop_hdr = fits.getheader(self._lyotstop_fname) result[0].header.set('MODE', self.mode, comment='Observing mode') result[0].header.set('CAMERA', self.camera, comment='Imager or IFS') result[0].header.set('APODIZER', self.apodizer, comment='Apodizer') result[0].header.set('APODTRAN', os.path.basename(self._apodizer_fname), comment='Apodizer transmission') result[0].header.set('PUPLSCAL', apodizer_hdr['PUPLSCAL'], comment='Apodizer pixel scale in m/pixel') result[0].header.set('PUPLDIAM', apodizer_hdr['PUPLDIAM'], comment='Full apodizer array size, incl padding.') result[0].header.set('FPM', self.fpm, comment='Focal plane mask') result[0].header.set('FPMTRAN', os.path.basename(self._fpm_fname), comment='FPM transmission') result[0].header.set('FPMSCAL', fpm_hdr['PIXSCALE'], comment='FPM spatial sampling, arcsec/pix') result[0].header.set('LYOTSTOP', self.lyotstop, comment='Lyot stop') result[0].header.set('LSTRAN', os.path.basename(self._lyotstop_fname), comment='Lyot stop transmission') result[0].header.set('PUPLSCAL', lyotstop_hdr['PUPLSCAL'], comment='Lyot stop pixel scale in m/pixel') result[0].header.set('PUPLDIAM', lyotstop_hdr['PUPLDIAM'], comment='Lyot stop array size, incl padding.')