In [1]:
# Makes print and division act like Python 3
from __future__ import print_function, division
# Import the usual libraries
import numpy as np
import matplotlib
import matplotlib.pyplot as plt
# Enable inline plotting at lower left
%matplotlib inline
from IPython.display import display, Latex, clear_output
In [2]:
# import main module
import pynrc
from pynrc import nrc_utils
# import pysynphot instance
#from pynrc.nrc_utils import S
Log messages for pynrc follow the same the logging functionality included in webbpsf. Logging levels include DEBUG, INFO, WARN, and ERROR.
In [3]:
pynrc.setup_logging()
If you get tired of the INFO level messages, simply type:
pynrc.setup_logging('WARN', verbose=False)
The basic NIRCam object consists of all the instrument settings one would specify for a JWST observation, including filter, pupil, and coronagraphic mask selections along with detector subarray settings and ramp sampling cadence (i.e., MULTIACCUM).
The NIRCam class makes use of high order polynomial coefficient maps to quickly generate large numbers of monochromatic PSFs that can be convolved with arbitrary spectra and collapsed into a final broadband PSF (or dispersed with NIRCam's slitless grisms). The PSF coefficients are calculated from a series of WebbPSF monochromatic PSFs and saved to disk. These polynomial coefficients are further modifed based on focal plane position and drift in the wavefront error relative to nominal OPD mao.
There are a multitude of posssible keywords one can pass upon initialization, including certain detector settings and PSF generation parameters. If not passed initially, then defaults are assumed. The user can update these parameters at any time by either setting attributes directly (e.g., filter, mask, pupil, etc.) along with using the update_detectors() and update_psf_coeff() methods.
For instance,
nrc = pynrc.NIRCam('F210M')
nrc.module = 'B'
nrc.update_detectors(read_mode='DEEP8', nint=10, ngroup=5)
is the same as:
nrc = pynrc.NIRCam('F210M', module='B', read_mode='DEEP8', nint=10, ngroup=5)
To start, we'll set up a simple observation using the F430M filter. Defaults will be populated for unspecified attributes such as module, pupil, mask, etc.
Check the function docstrings for more detailed information
In [4]:
nrc = pynrc.NIRCam(filter='F430M')
print('Filter: {}; Pupil: {}; Mask: {}; Module: {}'\
.format(nrc.filter, nrc.pupil, nrc.mask, nrc.module))
Keyword information for detector and PSF settings are stored in the det_info and psf_info dictionaries. These cannot be modified directly, but instead are updated via the update_detectors() and update_psf_coeff() methods.
In [5]:
print('Detector Info Keywords:')
print(nrc.det_info)
print('')
print('PSF Info Keywords:')
print(nrc.psf_info)
PSF coefficient information is stored in the psf_coeff attribute. This data is accessed by many of the NIRCam class functions to generate PSFs with arbitrary wavelength weights, such as the gen_psf() function.
In [8]:
# Demonstrate the color difference of the PSF for different spectral types, same magnitude
sp_M0V = pynrc.stellar_spectrum('M0V', 10, 'vegamag', nrc.bandpass)
sp_A0V = pynrc.stellar_spectrum('A0V', 10, 'vegamag', nrc.bandpass)
# Generate oversampled PSFs (counts/sec)
_, psf_M0V = nrc.gen_psf(sp_M0V, return_oversample=True)
_, psf_A0V = nrc.gen_psf(sp_A0V, return_oversample=True)
fig, axes = plt.subplots(1,3, figsize=(12,4))
axes[0].imshow(psf_M0V**0.5)
axes[0].set_title('M0V PSF ({})'.format(nrc.filter))
axes[1].imshow(psf_A0V**0.5)
axes[1].set_title('A0V PSF ({})'.format(nrc.filter))
diff = psf_M0V - psf_A0V
minmax = np.abs(diff).max() / 2
axes[2].imshow(diff, cmap='RdBu', vmin=-minmax, vmax=minmax)
axes[2].set_title('Difference')
fig.tight_layout()
Bandpass information is stored in the bandpass attribute and can be plotted with the convenience function plot_bandpass().
In [9]:
nrc.plot_bandpass()
Out[9]:
In [10]:
# Turn off those pesky informational texts
pynrc.setup_logging('WARN', verbose=False)
# Configure the observation for CDS frames (ngroup=2)
# Print out frame and ramp information using verbose=True
nrc.update_detectors(ngroup=2, verbose=True)
The sat_limits() function returns a dictionary of results. There's the option in include a Pysynphot spectrum, but if none is specificed the it defaults to a G2V star.
In [11]:
# Set verbose=True to print results in a user-friendly manner
sat_lims = nrc.sat_limits(verbose=True)
# Dictionary information
print("\nDictionary Info:", sat_lims)
By default, the function sat_limits() uses a G2V stellar spectrum, but any arbritrary spectrum can be passed via the sp keyword. In addition, using the bp_lim keyword, you can use spectral information to determine the brightness in some other bandpass that saturates the source within the NIRCam filter.
In [12]:
# Spectrum of an M0V star (not normalized)
sp_M0V = pynrc.stellar_spectrum('M0V')
# 2MASS Ks Bandpass
bp_k = pynrc.bp_2mass('K')
sat_lims = nrc.sat_limits(sp=sp_M0V, bp_lim=bp_k, verbose=True)
Now, let's get the same saturation limit assuming a 128x128 subarray.
In [13]:
nrc.update_detectors(wind_mode='WINDOW', xpix=128, ypix=128)
sat_lims = nrc.sat_limits(sp=sp_M0V, bp_lim=bp_k, verbose=True)
You can also use the saturation_levels() function to generate an image of a point source indicating the fractional well fill level.
In [14]:
# Spectum of A0V star with Ks = 8 mag
sp = pynrc.stellar_spectrum('M0V', 8, 'vegamag', bp_k)
sat_levels = nrc.saturation_levels(sp, full_size=False, ngroup=nrc.det_info['ngroup'])
print('Max Well Fraction: {:.2f}'.format(sat_levels.max()))
In [15]:
# Plot the well fill levels for each pixel
fig, ax = plt.subplots(1,1)
extent = 0.5*nrc.psf_info['fov_pix'] * np.array([-1,1,-1,1])
cax = ax.imshow(sat_levels, extent=extent, vmin=0, vmax=1)
ax.set_xlabel('Pixels')
ax.set_ylabel('Pixels')
ax.set_title('Well Fraction in {} of $K_s = 5$ M0V star'.format(nrc.filter))
cbar = fig.colorbar(cax)
cbar.set_label('Well Fill Fraction')
ax.tick_params(axis='both', color='white', which='both')
for k in ax.spines.keys():
ax.spines[k].set_color('white')
Information for slitless grism observations show wavelength-dependent results.
In [16]:
nrc = pynrc.NIRCam('F444W', pupil='GRISM0', ngroup=2, wind_mode='STRIPE', ypix=128)
sat_lims = nrc.sat_limits(sp=sp_M0V, bp_lim=bp_k, verbose=True)
Similarly, we can determine sensitivity limits of point sources (and extended sources) for the defined instrument configuration. By default, the sensitivity() function uses a flat spectrum. In this case, let's find the sensitivities NIRCam can reach in a single ~1000sec integration with the F430M filter. Noise values will depend on the exact MULTIACCUM settings.
In [17]:
nrc = pynrc.NIRCam('F430M')
nrc.update_detectors(read_mode='MEDIUM8', ngroup=10)
# The multiaccum_times attribute describes the various timing information
print(nrc.multiaccum_times)
In [18]:
sens = nrc.sensitivity(nsig=5, units='vegamag', verbose=True)
The sensitivity function also includes a keyword forwardSNR, which allows the user to pass a normalized spectrum and estimate the SNR For some extraction aperture.
In [19]:
sp = pynrc.stellar_spectrum('M0V', 20, 'vegamag', nrc.bandpass)
snr = nrc.sensitivity(sp=sp, forwardSNR=True, units='vegamag', verbose=True)
Armed with these two basic functions, we can attempt to determine the best instrument settings to optimize for SNR and efficiency. In these types of optimizations, we must consider observational constraints such as saturation levels, SNR requirements, and limits on acquisition time.
Note: The reported acquisition times do not include obsevatory and instrument-level overheads, such as slew times, filter changes, script compilations, etc. It only includes detector readout times (including reset frames and Fast Row Resets).
For instance, we want to observe an M-Dwarf (K=18 mag) in the F430M filter. What is the most efficient configuration to obtain an SNR of 100?
In [23]:
# Setup observation
nrc = pynrc.NIRCam('F430M', wind_mode='WINDOW', xpix=160, ypix=160)
# Spectrum of an M2V star
bp_k = pynrc.bp_2mass('K')
sp_M0V = pynrc.stellar_spectrum('M2V', 18, 'vegamag', bp_k)
In [24]:
# Run optimizer. Result is a ranked list sorted by efficiency.
tbl = nrc.ramp_optimize(sp_M0V, snr_goal=100, ng_min=5, nint_min=10, verbose=True)
For a slightly more complicated scenario, consider an additional foreground source. In this scenario, the F0V star will saturate much more quickly compared to the fainter M2V, so it limits which ramp settings we may want to use (assuming we want unsaturated frames, which isn't always necessarily true).
In [25]:
sp_F0V = pynrc.stellar_spectrum('F0V', 10, 'vegamag', bp_k)
tbl = nrc.ramp_optimize(sp_M0V, sp_bright=sp_F0V,
snr_goal=100, ng_min=5, nint_min=10, verbose=True)
If there are no objections to saturating the bright source, then we can set the well_frac_max parameter to something like 5 times the hard saturation limit. This allows for more efficient exposure settings.
In [26]:
tbl = nrc.ramp_optimize(sp_M0V, sp_bright=sp_F0V, well_frac_max=5,
snr_goal=100, ng_min=5, nint_min=10, verbose=True)
In [ ]: