Basic Usage

This tutorial walks through the basic usage of the pynrc package to calculate sensitivities and saturation limits for NIRCam in a variety of modes.

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

Getting Started

We assume you have already installed pynrc as outlined in the documentation.

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()
pyNRC log messages of level INFO and above will be shown.
pyNRC log outputs will be directed to the screen.

If you get tired of the INFO level messages, simply type:

pynrc.setup_logging('WARN', verbose=False)

First NIRCam Observation

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))
[     pynrc:INFO] Initializing SCA 485/A5
[     pynrc:INFO] RAPID readout mode selected.
[     pynrc:INFO] Setting nf=1, nd1=0, nd2=0, nd3=0.
[     pynrc:INFO] Updating PSF coeff with fov_pix=33 and oversample=4
Filter: F430M;  Pupil: CLEAR;  Mask: None;  Module: A

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)
Detector Info Keywords:
{'wind_mode': 'FULL', 'xpix': 2048, 'ypix': 2048, 'x0': 0, 'y0': 0, 'read_mode': 'RAPID', 'nint': 1, 'ngroup': 1, 'nf': 1, 'nd1': 0, 'nd2': 0, 'nd3': 0}

PSF Info Keywords:
{'fov_pix': 33, 'oversample': 4, 'offset_r': 0, 'offset_theta': 0, 'tel_pupil': None, 'save': True, 'force': False, 'opd': ('OPD_RevW_ote_for_NIRCam_requirements.fits', 0), 'jitter': 'gaussian', 'jitter_sigma': 0.007}

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 [6]:
# 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')
axes[1].imshow(psf_A0V**0.5)
axes[1].set_title('A0V PSF')

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()
../_images/tutorials_Basic_Usage_12_0.png

Bandpass information is stored in the bandpass attribute and can be plotted with the convenience function plot_bandpass().

In [7]:
nrc.plot_bandpass()
Out[7]:
<matplotlib.axes._subplots.AxesSubplot at 0x11c379208>
../_images/tutorials_Basic_Usage_14_1.png

1. Saturation Limits

One of the most basic functions is to determine the saturation limit of a CDS observation, so let’s try this for the current filter selection. Generally, saturation is considered to be 80% of the full well.

In [8]:
# 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)
New Ramp Settings:
  read_mode :    RAPID
  nf        :        1
  nd2       :        0
  ngroup    :        2
  nint      :        1
New Detector Settings
  wind_mode :     FULL
  xpix      :     2048
  ypix      :     2048
  x0        :        0
  y0        :        0
New Ramp Times
  t_group   :   10.737
  t_frame   :   10.737
  t_int     :   21.474
  t_int_tot :   32.210
  t_exp     :   21.474
  t_acq     :   32.216

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 [9]:
# Set verbose=True to print results in a user-friendly manner
sat_lims = nrc.sat_limits(verbose=True)

# Dictionary information
print(sat_lims)
F430M Saturation Limit assuming G2V source: 12.18 vegamag
{'satmag': 12.180657842243871, 'units': 'vegamag', 'Spectrum': 'G2V', 'bp_lim': 'F430M'}

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 [10]:
# 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)
Ks Saturation Limit for F430M assuming M0V source: 12.37 vegamag

Now, let’s get the same saturation limit assuming a 128x128 subarray.

In [11]:
nrc.update_detectors(wind_mode='WINDOW', xpix=128, ypix=128)
sat_lims = nrc.sat_limits(sp=sp_M0V, bp_lim=bp_k, verbose=True)
Ks Saturation Limit for F430M assuming M0V source: 7.94 vegamag

You can also use the saturation_levels() function to generate an image of a point source indicating the fractional well fill level.

In [12]:
# 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()))
Max Well Fraction: 0.76
In [13]:
# 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')
../_images/tutorials_Basic_Usage_25_0.png

Information for slitless grism observations show wavelength-dependent results.

In [14]:
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)
Ks Saturation Limit for F444W assuming M0V source:
   Wave   Sat Limit (vegamag)
--------- -------------------
     3.90                4.50
     4.00                4.49
     4.10                4.38
     4.20                4.26
     4.30                4.11
     4.40                3.87
     4.50                3.69
     4.60                3.53
     4.70                3.35
     4.80                3.18
     4.90                2.97
     5.00                2.38

2. Sensitivity Limits

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 [15]:
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)
{'t_frame': 10.73677, 't_group': 107.3677, 't_int': 1052.20346, 't_exp': 1052.20346, 't_acq': 1062.94547, 't_int_tot': 1062.94023}
In [16]:
sens = nrc.sensitivity(nsig=5, units='vegamag', verbose=True)
Point Source Sensitivity (5-sigma): 23.32 vegamag
Surface Brightness Sensitivity (5-sigma): 21.22 vegamag/arcsec^2

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 [17]:
sp = pynrc.stellar_spectrum('M0V', 20, 'vegamag', nrc.bandpass)
snr = nrc.sensitivity(sp=sp, forwardSNR=True, units='vegamag', verbose=True)
Point Source SNR (20.00 vegamag): 75.53 sigma
Surface Brightness SNR (20.00 vegamag/arcsec^2): 14.47 sigma

3. Ramp Optimization

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 [18]:
# 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 [19]:
# 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)
BRIGHT1
BRIGHT2
DEEP2
DEEP8
MEDIUM2
MEDIUM8
RAPID
SHALLOW2
SHALLOW4
 Pattern   NGRP NINT   t_int     t_exp     t_acq     SNR      Well     eff
---------- ---- ---- --------- --------- --------- -------- -------- --------
DEEP8         5   10     24.52    245.20    248.04    124.5    0.006    7.904
MEDIUM8       6   12     16.16    193.93    197.34    100.2    0.004    7.133
MEDIUM8       6   13     16.16    210.09    213.78    104.3    0.004    7.133
MEDIUM2       8   10     20.06    200.62    203.46     98.6    0.005    6.913
MEDIUM2       8   11     20.06    220.68    223.80    103.4    0.005    6.913
SHALLOW4     10   15     13.65    204.80    209.06     99.4    0.003    6.876
SHALLOW4     10   16     13.65    218.45    222.99    102.7    0.003    6.876
DEEP2         5   10     22.85    228.48    231.32    102.9    0.005    6.766
MEDIUM8       5   18     13.37    240.74    245.85    104.3    0.003    6.654
SHALLOW4      9   20     12.26    245.20    250.88    104.3    0.003    6.583
...         ...  ...       ...       ...       ...      ...      ...      ...
RAPID        10  511      2.79   1423.85   1568.87    101.1    0.001    2.551
RAPID        10  512      2.79   1426.64   1571.94    101.2    0.001    2.551
RAPID        10  513      2.79   1429.42   1575.01    101.3    0.001    2.551
RAPID        10  514      2.79   1432.21   1578.08    101.4    0.001    2.551
RAPID        10  515      2.79   1435.00   1581.15    101.5    0.001    2.551
RAPID        10  516      2.79   1437.78   1584.22    101.6    0.001    2.551
RAPID        10  517      2.79   1440.57   1587.29    101.7    0.001    2.551
RAPID        10  518      2.79   1443.36   1590.36    101.8    0.001    2.551
RAPID        10  519      2.79   1446.14   1593.43    101.8    0.001    2.551
RAPID        10  520      2.79   1448.93   1596.50    101.9    0.001    2.551
RAPID        10  521      2.79   1451.71   1599.57    102.0    0.001    2.551
Length = 50 rows

For a more complicated scenario, consider the M-Dwarf orbiting a F0V primary. In this scenario, the F0V star will saturate much more quickly compared to the fainter companion, so it limits which ramp settings we can use (assuming we want unsaturated frames, which isn’t always necessarily true).

In [20]:
sp_F0V = pynrc.stellar_spectrum('F0V', 13.7, 'vegamag', bp_k)
tbl = nrc.ramp_optimize(sp_M0V, sp_bright=sp_F0V,
                        snr_goal=100, ng_min=5, nint_min=10, verbose=True)
BRIGHT1
BRIGHT2
DEEP2
DEEP8
MEDIUM2
MEDIUM8
RAPID
SHALLOW2
SHALLOW4
 Pattern   NGRP NINT   t_int     t_exp     t_acq     SNR      Well     eff
---------- ---- ---- --------- --------- --------- -------- -------- --------
DEEP8         5   10     24.52    245.20    248.04    124.5    0.227    7.904
MEDIUM8       6   12     16.16    193.93    197.34    100.2    0.149    7.133
MEDIUM8       6   13     16.16    210.09    213.78    104.3    0.149    7.133
MEDIUM2       8   10     20.06    200.62    203.46     98.6    0.185    6.913
MEDIUM2       8   11     20.06    220.68    223.80    103.4    0.185    6.913
SHALLOW4     10   15     13.65    204.80    209.06     99.4    0.126    6.876
SHALLOW4     10   16     13.65    218.45    222.99    102.7    0.126    6.876
DEEP2         5   10     22.85    228.48    231.32    102.9    0.211    6.766
MEDIUM8       5   18     13.37    240.74    245.85    104.3    0.124    6.654
SHALLOW4      9   20     12.26    245.20    250.88    104.3    0.113    6.583
...         ...  ...       ...       ...       ...      ...      ...      ...
RAPID        10  511      2.79   1423.85   1568.87    101.1    0.026    2.551
RAPID        10  512      2.79   1426.64   1571.94    101.2    0.026    2.551
RAPID        10  513      2.79   1429.42   1575.01    101.3    0.026    2.551
RAPID        10  514      2.79   1432.21   1578.08    101.4    0.026    2.551
RAPID        10  515      2.79   1435.00   1581.15    101.5    0.026    2.551
RAPID        10  516      2.79   1437.78   1584.22    101.6    0.026    2.551
RAPID        10  517      2.79   1440.57   1587.29    101.7    0.026    2.551
RAPID        10  518      2.79   1443.36   1590.36    101.8    0.026    2.551
RAPID        10  519      2.79   1446.14   1593.43    101.8    0.026    2.551
RAPID        10  520      2.79   1448.93   1596.50    101.9    0.026    2.551
RAPID        10  521      2.79   1451.71   1599.57    102.0    0.026    2.551
Length = 50 rows