WFIRST test case fails if no pysynphot is present

[mperrin]

The WFIRST imager class fails for the fallback code for the no-pysynphot case. Looks like something minor in terms of the filter file setup or header, hopefully easy to fix.

I'm going to make it xfail for now, in the hopes of almost having the tests actually pass on Travis.

webbpsf/tests/test_wfirst.py:3: test_wfirstimager_psf FAILED

====================================================================================== FAILURES =======================================================================================
________________________________________________________________________________ test_wfirstimager_psf ________________________________________________________________________________

    def test_wfirstimager_psf():
        """Just test that instantiating WFIRSTImager works and can compute a PSF"""
        wi = wfirst.WFIRSTImager()
>       wi.calcPSF()

webbpsf/tests/test_wfirst.py:6:
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _

self = <webbpsf.wfirst.WFIRSTImager object at 0x1101dd610>, outfile = None, source = None, filter = None, nlambda = 10, monochromatic = None, fov_arcsec = 5.0, fov_pixels = None
oversample = 4, detector_oversample = 4, fft_oversample = 4, calc_oversample = None, rebin = True, clobber = True, display = False, return_intermediates = False, kwargs = {}
local_options = {'detector_oversample': 4, 'fft_oversample': 4, 'monochromatic': None, 'nlambda': 10}, fov_spec = 'arcsec = 5.000000'

    def calcPSF(self, outfile=None, source=None, filter=None, nlambda=None, monochromatic=None,
                fov_arcsec=None, fov_pixels=None,  oversample=None, detector_oversample=None,
                fft_oversample=None, calc_oversample=None, rebin=True, clobber=True, display=False,
                return_intermediates=False, **kwargs):
        """ Compute a PSF.

            The result can either be written to disk (set outfile="filename") or else will be returned as
            an astropy.io.fits HDUList object.


            Output sampling may be specified in one of two ways:

            1) Set `oversample=<number>`. This will use that oversampling factor beyond detector pixels
               for output images, and beyond Nyquist sampling for any FFTs to prior optical planes.
            2) set `detector_oversample=<number>` and `fft_oversample=<other_number>`. This syntax lets
               you specify distinct oversampling factors for intermediate and final planes. This is generally
               only relevant in the case of coronagraphic calculations.

            By default, both oversampling factors are set equal to 4. This default can be changed in your
            webbpsf configuration file.

            Notes
            -----
            More advanced PSF computation options (pupil shifts, source positions, jitter, ...)
            may be set by configuring the `.options` dictionary attribute of this class.

            Parameters
            ----------
            filter : string, optional
                Filter name. Setting this is just a shortcut for setting the object's filter first, then
                calling calcPSF afterwards.
            source : pysynphot.SourceSpectrum or dict or tuple
                specification of source input spectrum. Default is a 5700 K sunlike star.
            nlambda : int
                How many wavelengths to model for broadband?
                The default depends on how wide the filter is, as set by a lookup table in the webbpsf data distribution.
            monochromatic : float, optional
                Setting this to a wavelength value (in meters) will compute a monochromatic PSF at that
                wavelength, overriding filter and nlambda parameters.
            fov_arcsec : float
                field of view in arcsec. Default=5
            fov_pixels : int
                field of view in pixels. This is an alternative to fov_arcsec.
            outfile : string
                Filename to write. If None, then result is returned as an HDUList
            oversample, detector_oversample, fft_oversample : int
                How much to oversample. Default=4. By default the same factor is used for final output
                pixels and intermediate optical planes, but you may optionally use different factors
                if so desired.
            rebin : bool, optional
                If set, the output file will contain a FITS image extension containing the PSF rebinned
                onto the actual detector pixel scale. Thus, setting oversample=<N> and rebin=True is
                the proper way to obtain high-fidelity PSFs computed on the detector scale. Default is True.
            clobber : bool
                overwrite output FITS file if it already exists?
            display : bool
                Whether to display the PSF when done or not.
            save_intermediates, return_intermediates : bool
                Options for saving to disk or returning to the calling function the intermediate optical planes during the propagation.
                This is useful if you want to e.g. examine the intensity in the Lyot plane for a coronagraphic propagation. These have no
                effect for simple direct imaging calculations.


            For additional arguments, see the documentation for poppy.OpticalSystem.calcPSF()


            Returns
            -------
            outfits : fits.HDUList
                The output PSF is returned as a fits.HDUlist object.
                If `outfile` is set to a valid filename, the output is also written to that file.


            """

        if calc_oversample is not None:
            raise DeprecationWarning("The calc_oversample parameter is deprecated and will be removed in webbpsf 0.4. User fft_oversample instead.")
            fft_oversample = calc_oversample # back compatibility hook for deprecated arg name.

        _log.info("Setting up PSF calculation for "+self.name)

        # first make sure that webbpsf's configuration is used to override any of the
        # same configuration options in poppy. This is admittedly perhaps overbuilt to have identical
        # settings in both packages, but the intent is to shield typical users of webbpsf
        # from having to think about the existence of the underlying library. They can
        # just deal with one set of settings.
        #config._apply_settings_to_poppy()

        if filter is not None:
            self.filter = filter

        local_options = self.options.copy()  # all local state should be stored in a dict, for
                                      # ease of handing off to the various subroutines of
                                      # calcPSF. Don't just modify the global self.options
                                      # structure since that would pollute it with temporary
                                      # state as well as persistent state.
        local_options['monochromatic'] = monochromatic



        #----- choose # of wavelengths intelligently. Do this first before generating the source spectrum weighting.
        if nlambda is None or nlambda==0:
            # Automatically determine number of appropriate wavelengths.
            # Make selection based on filter configuration file
            try:
                nlambda = self._filter_nlambda_default[self.filter]
                _log.debug("Automatically selecting # of wavelengths: %d" % nlambda)
            except:
                nlambda=10
                _log.warn("Filter %s not found in lookup table for default number of wavelengths to use.. setting default nlambda=%d" % (self.filter, nlambda))
        local_options['nlambda'] = nlambda



        #----- calculate field of view depending on supplied parameters
        if fov_arcsec is None and fov_pixels is None:  #pick decent defaults.
            if self.name =='MIRI': fov_arcsec=12.
            else: fov_arcsec=5.
            fov_spec = 'arcsec = %f' % fov_arcsec
        elif fov_pixels is not None:

            if np.isscalar(fov_pixels):
                fov_spec = 'pixels = %d' % fov_pixels
            else:
                fov_spec = 'pixels = (%d, %d)' % (fov_pixels[0], fov_pixels[1])
        elif fov_arcsec is not None:
            if np.isscalar(fov_arcsec):
                fov_spec = 'arcsec = %f' % fov_arcsec
            else:
                fov_spec = 'arcsec = (%.3f, %.3f)' % (fov_arcsec[0], fov_arcsec[1])

        _log.debug('FOV set to '+fov_spec)

        #---- Implement the semi-convoluted logic for the oversampling options. See docstring above
        if oversample is not None and detector_oversample is not None and fft_oversample is not None:
            # all options set, contradictorily -> complain!
            raise ValueError("You cannot specify simultaneously the oversample= option with the detector_oversample and fft_oversample options. Pick one or the other!")
        elif oversample is None and detector_oversample is None and fft_oversample is None:
            # nothing set -> set oversample = 4
            oversample = conf.default_oversampling
        if detector_oversample is None: detector_oversample = oversample
        if fft_oversample is None: fft_oversample = oversample
        local_options['detector_oversample']=detector_oversample
        local_options['fft_oversample']=fft_oversample

        #----- compute weights for each wavelength based on source spectrum
>       wavelens, weights = self._getWeights(source=source, nlambda=nlambda, monochromatic=monochromatic)

webbpsf/webbpsf_core.py:398:
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _

self = <webbpsf.wfirst.WFIRSTImager object at 0x1101dd610>, source = None, nlambda = 10, monochromatic = None, verbose = False

    def _getWeights(self, source=None, nlambda=5, monochromatic=None, verbose=False):
        """ Return the set of discrete wavelengths, and weights for each wavelength,
            that should be used for a PSF calculation.

            Uses pysynphot (if installed), otherwise assumes simple-minded flat spectrum

            """
        if monochromatic is not None:
            poppy_core._log.info(" monochromatic calculation requested.")
            return (np.asarray([monochromatic]),  np.asarray([1]) )

        elif _HAS_PYSYNPHOT and (isinstance(source, pysynphot.spectrum.SourceSpectrum)  or source is None):
            """ Given a pysynphot.SourceSpectrum object, perform synthetic photometry for
                nlambda bins spanning the wavelength range of interest.

                Because this calculation is kind of slow, cache results for reuse in the frequent
                case where one is computing many PSFs for the same spectral source.
                """
            poppy_core._log.debug("Calculating spectral weights using pysynphot, nlambda=%d, source=%s" % (nlambda, str(source)))
            if source is None:
                try:
                    source = pysynphot.Icat('ck04models',5700,0.0,2.0)
                except:
                    poppy_core._log.error("Could not load Castelli & Kurucz stellar model from disk; falling back to 5700 K blackbody")
                    source = pysynphot.BlackBody(5700)
            poppy_core._log.debug("Computing spectral weights for source = "+str(source))

            try:
                key = self._getSpecCacheKey(source, nlambda)
                if key in self._spectra_cache.keys():
                    poppy_core._log.debug("Previously computed spectral weights found in cache, just reusing those")
                    return self._spectra_cache[keys]
            except:
                pass  # in case sourcespectrum lacks a name element so the above lookup fails - just do the below calc.

            poppy_core._log.info("Computing wavelength weights using synthetic photometry for %s..." % self.filter)
            band = self._getSynphotBandpass(self.filter)
            # choose reasonable min and max wavelengths
            w_above10 = np.where(band.throughput > 0.10*band.throughput.max())

            minwave = band.wave[w_above10].min()
            maxwave = band.wave[w_above10].max()
            poppy_core._log.debug("Min, max wavelengths = %f, %f" % (minwave/1e4, maxwave/1e4))

            wave_bin_edges =  np.linspace(minwave,maxwave,nlambda+1)
            wavesteps = (wave_bin_edges[:-1] +  wave_bin_edges[1:])/2
            deltawave = wave_bin_edges[1]-wave_bin_edges[0]
            effstims = []

            for wave in wavesteps:
                poppy_core._log.debug("Integrating across band centered at %.2f microns with width %.2f" % (wave/1e4,deltawave/1e4))
                box = pysynphot.Box(wave, deltawave) * band
                if box.throughput.max() == 0:  # watch out for pathological cases with no overlap (happens with MIRI FND at high nlambda)
                    result = 0.0
                else:
                    binset =  np.linspace(wave-deltawave, wave+deltawave, 30)  # what wavelens to use when integrating across the sub-band?
                    result = pysynphot.Observation(source, box, binset=binset).effstim('counts')
                effstims.append(result)

            effstims = np.array(effstims)
            effstims /= effstims.sum()
            wave_m =  band.waveunits.Convert(wavesteps,'m') # convert to meters

            newsource = (wave_m, effstims)
            if verbose: _log.info( " Wavelengths and weights computed from pysynphot: "+str( newsource))
            self._spectra_cache[ self._getSpecCacheKey(source,nlambda)] = newsource
            return newsource
        elif isinstance(source, dict) and ('wavelengths' in source) and ('weights' in source):
            # Allow providing directly a set of specific weights and wavelengths, as in poppy.calcPSF source option #2
            return (source['wavelengths'], source['weights'])
        elif isinstance(source, tuple) and len(source) == 2:
            # Allow user to provide directly a tuple, as in poppy.calcPSF source option #3
            return source

        else:  #Fallback simple code for if we don't have pysynphot.
            poppy_core._log.warning("Pysynphot unavailable (or invalid source supplied)!   Assuming flat # of counts versus wavelength.")
            # compute a source spectrum weighted by the desired filter curves.
            # The existing FITS files all have wavelength in ANGSTROMS since that is the pysynphot convention...
            filterfits = fits.open(self._filters[self.filter].filename)
            filterdata = filterfits[1].data
            try:
                f1 = filterdata.WAVELENGTH
                d2 = filterdata.THROUGHPUT
            except:
                raise ValueError("The supplied file, %s, does not appear to be a FITS table with WAVELENGTH and THROUGHPUT columns." % self._filters[self.filter].filename )
            if 'WAVEUNIT' in  filterfits[1].header.keys():
                waveunit  = filterfits[1].header['WAVEUNIT']
            else:
                poppy_core._log.warn("CAUTION: no WAVEUNIT keyword found in filter file {0}. Assuming = Angstroms by default".format(filterfits.filename()))
                waveunit = 'Angstrom'
>           if waveunit != 'Angstrom': raise ValueError("The supplied file, %s, does not have WAVEUNIT = Angstrom as expected." % self._filters[self.filter].filename )
E           ValueError: The supplied file, /Users/mperrin/software/webbpsf-data/WFIRSTImager/filters/Z087_throughput.fits, does not have WAVEUNIT = Angstrom as expected.

/Users/mperrin/Dropbox (Personal)/Documents/software/poppy/poppy/instrument.py:722: ValueError
----------------------------------------------------------------------------------- Captured stderr -----------------------------------------------------------------------------------
poppy     : INFO     Optic from /Users/mperrin/software/webbpsf-data/AFTA_symmetrical.fits: Loaded amplitude transmission from /Users/mperrin/software/webbpsf-data/AFTA_symmetrical.fits
webbpsf   : INFO     Setting up PSF calculation for WFIRSTImager
webbpsf   : WARNING  Filter Z087 not found in lookup table for default number of wavelengths to use.. setting default nlambda=10
poppy     : WARNING  Pysynphot unavailable (or invalid source supplied)!   Assuming flat # of counts versus wavelength.
================================================================== 1 failed, 36 passed, 1 xfailed in 361.24 seconds ==================================================================

[josePhoenix]

Ah, this is just a case of making the comparison case insensitive. I'll open a pull request to POPPY.

[josePhoenix]

Fixed through changes in POPPY.