In [1]:
import poppy
poppy.conf.autosave_fftw_wisdom = False  # workaround for recent versions of POPPY and Python 3.5+
import webbpsf
import webbpsf.wfirst
import numpy as np

webbpsf.setup_logging()
wfi = webbpsf.wfirst.WFI()

# wavelength (to convert coefficients into meters OPD)
wavelength = 0.5e-6

# coefficients in waves of OPD
coefficients_waves = np.array([ 
    0.00000000e+00,   0.00000000e+00,   0.00000000e+00, # Z1-3 (zero/omitted)
   -0.0469, -0.0037, 0.0105 # Z4-6
])
# coefficients in meters, as needed by poppy.ZernikeWFE
coefficients_meters = wavelength * coefficients_waves

zernike_optic = poppy.ZernikeWFE(radius=wfi.PUPIL_RADIUS, coefficients=coefficients_meters)


import types

def _get_aberrations(self):
   return zernike_optic  # defined in cell above

_get_aberrations = types.MethodType(_get_aberrations, wfi)

wfi._get_aberrations = _get_aberrations

# Calculate the PSF with the new set of aberrations
psf = wfi.calcPSF(monochromatic=0.5e-6, oversample = 44)


/Users/rubind/anaconda/envs/astroconda/lib/python2.7/site-packages/pysynphot/locations.py:119: UserWarning: Extinction files should be moved to $PYSYN_CDBS/extinction for compatibility with future versions of pysynphot.
  warnings.warn('Extinction files should be moved to '
/Users/rubind/anaconda/envs/astroconda/lib/python2.7/site-packages/pysynphot/locations.py:155: UserWarning: Extinction files not found in /Users/rubind/Downloads/synphot1-2/grid/extinction
  warnings.warn('Extinction files not found in %s' % (extdir,))
WebbPSF log messages of level INFO and above will be shown.
WebbPSF log outputs will be directed to the screen.
[webbpsf] Setting up PSF calculation for WFI
[  poppy] Monochromatic calculation requested.
[webbpsf] Using the unmasked WFI pupil shape based on wavelengths requested
[webbpsf] PSF calc using fov_arcsec = 5.000000, oversample = 44, number of wavelengths = 1
[webbpsf] Creating optical system model:
[  poppy] Initialized OpticalSystem: WFIRST+WFI
[  poppy] WFIRST Entrance Pupil: Loaded amplitude transmission from /Users/rubind/anaconda/envs/astroconda/share/webbpsf-data/wfc_pupil_rev_mcr.fits
[  poppy] WFIRST Entrance Pupil: Loaded OPD from /Users/rubind/anaconda/envs/astroconda/share/webbpsf-data/upscaled_HST_OPD.fits
[  poppy] Added pupil plane: WFIRST Entrance Pupil
[  poppy] Added coordinate inversion plane: OTE exit pupil
[  poppy] Added pupil plane: Zernike WFE
[  poppy] Added detector with pixelscale=0.11 and oversampling=44: WFI detector
[  poppy] Calculating PSF with 1 wavelengths
[  poppy]  Propagating wavelength = 5e-07 m
[  poppy]   Calculation completed in 4.897 s
[  poppy] PSF Calculation completed.
[  poppy]  Adding extension with image downsampled to detector pixel scale.

In [2]:
psf.info()


Filename: (No file associated with this HDUList)
No.    Name         Type      Cards   Dimensions   Format
0    OVERSAMP    PrimaryHDU      42   (1980, 1980)   float64   
1    DET_SAMP    ImageHDU        44   (45, 45)     float64   

In [3]:
psf.writeto("IFC_PSF_0.5_micron_2.5_mas.fits", clobber=True)


WARNING: VerifyWarning: Card is too long, comment will be truncated. [astropy.io.fits.card]
[astropy] VerifyWarning: Card is too long, comment will be truncated.

In [1]:
import pyfits
f = pyfits.open("IFC_PSF_0.5_micron_2.5_mas.fits")
newdat = f[0].data[0::2,0::2] + f[0].data[1::2,0::2] + f[0].data[0::2,1::2] + f[0].data[1::2,1::2]
f[0].data = newdat
f.writeto("IFC_PSF_0.5_micron_5_mas.fits", clobber=True)
f.close()
newdat.shape


/Users/rubind/anaconda/envs/astroconda/lib/python2.7/site-packages/pyfits/file.py:372: UserWarning: Overwriting existing file 'IFC_PSF_0.5_micron_5_mas.fits'.
  warnings.warn("Overwriting existing file %r." % self.name)
Out[1]:
(990, 990)

In [ ]:
from scipy.interpolate import RectBivariateSpline
from scipy import fftpack as ft
import pyfits
import numpy as np


dat = np.zeros([41, 41], dtype=np.float64)
dat[20,20]=1
vals = np.arange(41.) - 20

ifn = RectBivariateSpline(vals, vals, dat, kx = 3, ky = 3, s = 0)
spl_space = 0.05

for rotang in [0, np.pi/8., np.pi/4., np.pi/2., np.pi*(2./3.), np.pi*(4./3.)]:
    newvals = np.arange(990.)*(spl_space/0.05)/10.
    newvals -= np.median(newvals)
    xvals = newvals
    yvals = newvals


    basis = np.zeros([len(xvals), len(yvals)], dtype=np.float64)
    for i in range(len(xvals)):
        if i % 100 == 0:
            print i
        for j in range(len(yvals)):
            xp = np.cos(rotang)*xvals[i] - np.sin(rotang)*yvals[j]
            yp = np.sin(rotang)*xvals[i] + np.cos(rotang)*yvals[j]
        
            basis[i,j] = ifn(xp, yp)
        
    basis.shape
    f = pyfits.open("IFC_PSF_0.5_micron_5_mas.fits")
    f[0].data = basis
    f.writeto("spline_rot=%.3f_spl=%.3f.fits" % (rotang, spl_space), clobber=True)
    adjust = basis*0
    for i in range(10):
        adjust[len(adjust)/2 + i - 5, len(adjust)/2] = 1.
    conv = np.real(ft.ifft2(ft.fft2(basis) * ft.fft2(newdat) * ft.fft2(adjust)))
    f[0].data = conv
    f.writeto("basis_5_mas_rot=%.3f_spl=%.3f.fits" % (rotang, spl_space), clobber=True)

    f.close()
print "Done!"


0
100
200
300
400
500
600
700
800
900
/Users/rubind/anaconda/envs/astroconda/lib/python2.7/site-packages/pyfits/file.py:372: UserWarning: Overwriting existing file 'spline_rot=0.000_spl=0.050.fits'.
  warnings.warn("Overwriting existing file %r." % self.name)
/Users/rubind/anaconda/envs/astroconda/lib/python2.7/site-packages/pyfits/file.py:372: UserWarning: Overwriting existing file 'basis_5_mas_rot=0.000_spl=0.050.fits'.
  warnings.warn("Overwriting existing file %r." % self.name)
0
100
200
300
400
500
600

In [ ]:
import pyfits
from scipy.interpolate import RectBivariateSpline
from scipy import fftpack as ft
import pyfits
import numpy as np



f = pyfits.open("IFC_PSF_0.5_micron_2.5_mas.fits")
newdat = f[0].data[0::2,0::2] + f[0].data[1::2,0::2] + f[0].data[0::2,1::2] + f[0].data[1::2,1::2]
f[0].data = newdat
f.writeto("IFC_PSF_0.5_micron_5_mas.fits", clobber=True)
f.close()
newdat.shape





dat = np.zeros([41, 41], dtype=np.float64)
dat[20,20]=1
vals = np.arange(41.) - 20

ifn = RectBivariateSpline(vals, vals, dat, kx = 3, ky = 3, s = 0)
spl_space = 0.05

for rotang in [0, np.pi/8., np.pi/4., np.pi/2., np.pi*(2./3.), np.pi*(4./3.)]:
    newvals = np.arange(990.)*(spl_space/0.05)/10.
    newvals -= np.median(newvals)
    xvals = newvals
    yvals = newvals


    basis = np.zeros([len(xvals), len(yvals)], dtype=np.float64)
    for i in range(len(xvals)):
        if i % 100 == 0:
            print i
        for j in range(len(yvals)):
            xp = np.cos(rotang)*xvals[i] - np.sin(rotang)*yvals[j]
            yp = np.sin(rotang)*xvals[i] + np.cos(rotang)*yvals[j]
        
            basis[i,j] = ifn(xp, yp)
        
    basis.shape
    f = pyfits.open("IFC_PSF_0.5_micron_5_mas.fits")
    f[0].data = basis
    f.writeto("spline_rot=%.3f_spl=%.3f.fits" % (rotang, spl_space), clobber=True)
    adjust = basis*0
    for i in range(10):
        adjust[len(adjust)/2 + i - 5, len(adjust)/2] = 1.
    conv = np.real(ft.ifft2(ft.fft2(basis) * ft.fft2(newdat) * ft.fft2(adjust)))
    f[0].data = conv
    f.writeto("basis_5_mas_rot=%.3f_spl=%.3f.fits" % (rotang, spl_space), clobber=True)

    f.close()
print "Done!"

In [ ]: