STIPS Tutorials
The Space Telescope Imaging Product Simulator (STIPS) is a versatile tool designed to simulate exposure level astronomical scenes from the Wide Field Instrument (WFI) for Roman. This article presents a series of tutorials demonstrating the range of functionalities, allowing users to develop simulations tailored to their scientific goals.
Tutorial #1: Basic STIPS Usage
The following tutorial describes a straightforward use of the STIPS software that can serve as the basis for other use cases. The Basic STIPS Usage Tutorial builds on the concepts introduced in the STIPS Overview article and is designed to walk through the phases of using STIPS at the most introductory level: creating a small scene, designing an observation, and generating a simulated image.
A. Generating a Simple Astronomical Scene
STIPS provides functionalities to generate catalogs of stars or galaxies with user-specified input parameters. Users can also provide their own source catalogs, as described in Catalogs formatting on readthedocs.
In the example below, we show how to create an astronomical scene by passing the user-defined input catalogs to STIPS in the FITS format. The catalog contains the coordinates of the sources, the observed flux, and any necessary shape parameters. The code block below provides an example of the required input information to generate a catalog containing two point sources. The catalog is saved to a file called catalog.fits for later use.
from astropy.io import fits
cols=[]
cols.append(fits.Column(name='id' ,array=[1,2] , format='K' )) # Object ID
cols.append(fits.Column(name='ra' ,array=[90.02,90.03] , format='D' )) # RA in degrees
cols.append(fits.Column(name='dec' ,array=[29.98,29.97] , format='D' )) # DEC in degrees
cols.append(fits.Column(name='flux' ,array=[0.00023,0.0004] , format='D' )) # Flux in `units`
cols.append(fits.Column(name='type' ,array=['point','point'], format='8A')) # `point` or `sersic`
cols.append(fits.Column(name='n' ,array=[0,0] , format='D' )) # Sersic profile index
cols.append(fits.Column(name='re' ,array=[0,0] , format='D' )) # Half-light radius in pixels
cols.append(fits.Column(name='phi' ,array=[0,0] , format='D' )) # Angle of PA in degrees
cols.append(fits.Column(name='ratio',array=[0,0] , format='D' )) # Axial Ratio
cols.append(fits.Column(name='notes',array=['',''] , format='8A')) # Notes
cols.append(fits.Column(name='units',array=['j','j'] , format='8A')) # Units, 'j' for jansky
# Create output fits table
hdut = fits.BinTableHDU.from_columns(cols)
hdut.header['TYPE']='mixed'
hdut.header['FILTER']='F129'
# Write to disk
hdut.writeto('catalog.fits',overwrite=True)
B. Initial Observation Setup
Next, we prepare inputs for the scene we created above. We use the ObservationModule to generate an observation object. We generate an obs dictionary containing the properties of the observation (e.g. the instrument and the detector setups) and parameters describing central coordinates (RA, Dec), orientation angle, and computational parameters for the simulation.
The Observation Dictionary
An example obs dictionary is specified in the code block below. In this example, we will use
- The imaging filter F129
- The WFI01 Sensor Chip Assembly (SCA)
- Sky background from Pandeia
- An observation ID of 42
- An exposure time of 300 seconds
Together with the following single offset setting:
- An offset ID of 1
- No centering
- No changes in RA, Dec, and PA from the center of the observation
from stips.observation_module import ObservationModule
# Build observation parameters
obs = {'instrument' : 'WFI',
'filters' : ['F129'],
'detectors' : 1,
'background' : 'pandeia',
'observations_id' : 42,
'exptime' : 300,
'offsets' : [{'offset_id' : 1 ,
'offset_centre': False,
'offset_ra' : 0.0 ,
'offset_dec' : 0.0 ,
'offset_pa' : 0.0 }]}
The Observation Object
An observation object combines the obs dictionary with the central coordinates (RA, Dec), orientation angle, and computational parameters for the simulation. Then, the setup can be initialized.
# Create observation object
obm = ObservationModule(obs,
ra = 90,
dec = 30,
pa = 0,
seed = 42,
cores = 6)
C. Generating a Simulated Image
The code block below combines the inputs required for
STIPS
that were generated above and finalizes the
ObservationModule
to run a simulation. The final simulated image will be saved under the name specified by the
fits_file
variable, and the result is shown on the right.
import numpy as np
import matplotlib.pyplot as plt
# Initialize the local instrument
obm.nextObservation()
# Add catalog with sources
cat_name = obm.addCatalogue('catalog.fits')
# Add error to image
psf_file = obm.addError()
# Call the final method
fits_file, mosaic_file, params = obm.finalize(mosaic=False)
print("Output FITS file is {}".format(fits_file))
data = fits.open(fits_file)[1].data
# Plot the simulated FITS image
fig, ax = plt.subplots(2, 1, figsize=(16, 10))
vmin = np.percentile(data, 5)
vmax = np.percentile(data, 95)
ax[0].imshow(data,
origin="lower",
cmap="grey",
vmin=vmin,
vmax=vmax)
ax[1].imshow(data,
origin="lower",
cmap="grey",
vmin=vmin,
vmax=vmax)
ax[1].set_xlim(2300, 3300)
ax[1].set_ylim(700, 1700)
ax[0].set_xlabel("X [pixels]")
ax[0].set_ylabel("Y [pixels]")
ax[1].set_xlabel("X [pixels]")
ax[1].set_ylabel("Y [pixels]")
plt.suptitle("Output Simulated FITS image")
ax[1].set_title("Zoomed portion of simulated image")
plt.tight_layout()
plt.show()
Result of STIPS Basic Use Tutorial
This output FITS image was simulated with STIPS using this tutorial.
Tutorial #2: Generating a Scene from a Simulated Catalog
In this example, we show users how to build a simulated catalog using STIPS built-in modules. The user must set the parameters for the stellar and galactic populations in the scene.
A. Generating an Astronomical Scene
Stellar Source Catalog
Below is an example showing how to generate a stellar population defined by stellar_parameters.
obs_prefix = 'notebook_example'
obs_ra = 150.0
obs_dec = -2.5
from stips.scene_module import SceneModule
scm = SceneModule(out_prefix=obs_prefix, ra=obs_ra, dec=obs_dec)
stellar_parameters = {
'n_stars': 100,
'age_low': 7.5e12,
'age_high': 7.5e12,
'z_low': -2.0,
'z_high': -2.0,
'imf': 'salpeter',
'alpha': -2.35,
'binary_fraction': 0.1,
'clustered': True,
'distribution': 'invpow',
'radius': 100.0,
'radius_units': 'pc',
'distance_low': 20.0,
'distance_high': 20.0,
'offset_ra': 0.0,
'offset_dec': 0.0
}
stellar_cat_file = scm.CreatePopulation(stellar_parameters)
print("Stellar population saved to file {}".format(stellar_cat_file))
Galaxy Source Catalog
Below is an example showing how to generate a galaxy population defined by galaxy_parameters.
galaxy_parameters = {
'n_gals': 10,
'z_low': 0.0,
'z_high': 0.2,
'rad_low': 0.01,
'rad_high': 2.0,
'sb_v_low': 30.0,
'sb_v_high': 25.0,
'distribution': 'uniform',
'clustered': False,
'radius': 200.0,
'radius_units': 'arcsec',
'offset_ra': 0.0,
'offset_dec': 0.0,
}
galaxy_cat_file = scm.CreateGalaxies(galaxy_parameters)
print("Galaxy population saved to file {}".format(galaxy_cat_file))
B. Generating the Scene
The stellar source catalog can now be observed. In this example, we will use:
- The F129 filter
- The WFI03 SCA
- Distortion on
- A uniform sky background level of 0.24 counts/s/pixel
- An ID of 1
- An exposure time of 1500 seconds
Together with the following single offset setting:
- An ID of 1
- No centering
- Changes in RA of 2 degrees, and PA of 0.5 degrees from the center of the observation
- No change in DEC
And the following residual settings:
- Flatfield residual on
- Dark current residual on
- Cosmic ray removal residual off
- Poisson noise residual off
- Readnoise residual off
# Set the parameters for the offset, residuals, and observation
offset_1 = {
'offset_id': 1,
'offset_centre': False,
'offset_ra': 2.0,
'offset_dec': 0.0,
'offset_pa': 0.5
}
residuals_1 = {
'residual_flat': True,
'residual_dark': True,
'residual_cosmic': False,
'residual_poisson': False,
'residual_readnoise': False
}
observation_parameters_1 = {
'instrument': 'WFI',
'filters': ['F129'],
'detectors': 3,
'distortion': True,
'background': 0.24,
'observations_id': 1,
'exptime': 1500,
'offsets': [offset_1]
}
# Create observation object
obm_1 = ObservationModule(observation_parameters_1,
out_prefix=obs_prefix,
ra=obs_ra,
dec=obs_dec,
residual=residuals_1)
# nextObservation is called to move between different combinations of offset and filter.
# It must be called once in order to initialize the observation module to the first observation before adding catalogs.
obm_1.nextObservation()
# Add the scene catalogues
output_stellar_catalogues = obm_1.addCatalogue(stellar_cat_file)
output_galaxy_catalogues = obm_1.addCatalogue(galaxy_cat_file)
print("Output Catalogues are {} and {}".format(output_stellar_catalogues, output_galaxy_catalogues))
# Add in the error residuals
obm_1.addError()
# Finalize the observation
fits_file, mosaic_file, params = obm_1.finalize(mosaic=False)
print("Output FITS file is {}".format(fits_file))
print("Output Mosaic file is {}".format(mosaic_file))
print("Observation Parameters are {}".format(params))
data = fits.open(fits_file)[1].data
vmin = np.percentile(data, 5)
vmax = np.percentile(data, 95)
plt.imshow(data,
origin="lower",
vmin=vmin,
vmax=vmax)
plt.xlabel("X [pixels]")
plt.ylabel("Y [pixels]")
plt.show()
Result of the Astronomical Scene Tutorial
Result of Tutorial #2. Image may look different due to random seed values on the source generators.
Tutorial #3: Modifying Observation Parameters
We can use the results of Tutorial #2 as the first of a set of observations with different properties (e.g., dither, noise residuals).
The function nextObservation()is required to move between different observations. This function must be called at least once, in order to initialize the observation module to the first observation before catalogs are added. The
obs_prefix
,
obs_ra
, and
obs_dec
values are set during initial setup of
STIPS
. Then, we can modify the offsets so that the same underlying astronomical scene is observed with different parameters. The code below adds an offset of 10 degrees in the RA and a rotation of 27 degrees, with a read noise residual instead of a dark current residual. This process can be repeated to generate sets of observations with different observing parameters, for example, to construct complex mosaics or compare results with different noise residuals.
# Set the parameters for the offset, residuals, and observation
offset_2 = {
'offset_id': 2,
'offset_centre': False,
'offset_ra': 10.0,
'offset_dec': 0.0,
'offset_pa': 27
}
residuals_2 = {
'residual_flat': True,
'residual_dark': False,
'residual_cosmic': False,
'residual_poisson': False,
'residual_readnoise': True
}
observation_parameters_2 = {
'instrument': 'WFI',
'filters': ['F129'],
'detectors': 3,
'distortion': True,
'background': 0.24,
'observations_id': 1,
'exptime': 1500,
'offsets': [offset_2]
}
# Create observation object
obm_2 = ObservationModule(observation_parameters_2,
out_prefix=obs_prefix,
ra=obs_ra, dec=obs_dec,
residual=residuals_2)
# Initialize the local instrument
obm_2.nextObservation()
Tutorial #4: PSFs and Adding Sources
STIPS
offers a variety of customizable options through the use of utility functions. As an example, we describe the use of the makePSF module. This module transforms the input PSFs from
WebbPSF
into effective PSFs, which are empirical models defining what fraction of a point source's flux will land in a particular pixel. The generated PSF models can then be used to create new scenes or to add sources to existing scenes.
A. Generating effective PSF models
STIPS
enables users to manually produce PSF models for a given set of detector positions and use them to insert point sources into an existing scene. PSFs at specific locations can be generated by providing the input PSF library and detector positions to the makePSF module within
STIPS
. The generated PSFs are returned as a 2D array.
STIPS performs two types of interpolation to obtain the exact PSF model at user-specified locations. First, it applies bilinear interpolation of the 3x3 PSF library array to compute the best PSF at the specified integer SCA pixels. Then, it performs bicubic interpolations over the supersampled pixel grid of this PSF to generate the best PSF at the specified SCA sub-pixel positions.
As an example, the file psf_WFI_2.0.0_F129_sca01.fits contains the library PSF for filter F129 and SCA 1. This file is available in the notebooks directory of the STIPS GitHub repository.
The library PSF can be loaded using the following code:
import stips
with fits.open('psf_WFI_2.0.0_F129_sca01.fits') as hdul:
test_psf = stips.utilities.makePSF.make_epsf(hdul[0].data[0])
vmin = np.percentile(test_psf, 5)
vmax = np.percentile(test_psf, 95)
plt.imshow(test_psf,
origin='lower',
vmin=vmin,
vmax=vmax)
plt.xlabel('X [pixels]')
plt.ylabel('Y [pixels]')
plt.show()
Visualization of the Example PSF
Visualization of one of the library PSF from the FITS file psf_WFI_2.0.0_F129_sca01.fits .
B. Adding an Artificial Point Source
We will now insert an artificial star into the astronomical scene defined in Example #2. For convenience, we provide the code used in Example #2 to generate the astronomical scene in the expandable box below.
The library PSF can be added to the simulated scene by specifying the source addition parameters, Location (px) and Flux in Jansky, as demonstrated in the code below. The result of this code applied to the scene from Example #2 is shown on the right.
import matplotlib.patches as patches
# Image created in the example above
with fits.open('notebook_example_1_0.fits') as result_file:
result_data = result_file[1].data
# Open library PSF generated earlier in this section
with fits.open('psf_WFI_2.0.0_F129_sca01.fits') as hdul:
test_psf = stips.utilities.makePSF.make_epsf(hdul[0].data[0])
# Default upscaling factor for STIPS
PSF_UPSCALE = 4
# Calculate center and size
psf_middle = int(np.round((test_psf[0].shape[0]-1) / 2))
boxsize = int(np.floor(psf_middle/PSF_UPSCALE))
added_source = stips.utilities.makePSF.place_source(2000,
2000,
3000,
result_data,
test_psf,
boxsize=boxsize,
psf_center=psf_middle)
# We will create a plot to compare the scene before+after adding the PSF
fig, ax = plt.subplots(2, 1, figsize=(8, 10))
ax[0].imshow(fits.getdata('notebook_example_1_0.fits'),
vmax=1,
origin="lower")
ax[1].imshow(added_source,
vmax=1,
origin="lower")
circle1 = patches.Circle((2000, 2000), 75,
fill=False,
edgecolor='red',
linewidth=1)
circle2 = patches.Circle((2000, 2000), 75,
fill=False,
edgecolor='red',
linewidth=1)
ax[0].add_patch(circle1)
ax[1].add_patch(circle2)
for subplot in ax.ravel():
subplot.set_xlabel("X [pixels]")
subplot.set_ylabel("Y [pixels]")
subplot.set_xlim(1500, 2500)
subplot.set_ylim(1500, 2500)
ax[0].set_title("Before adding PSF")
ax[1].set_title("After adding PSF")
plt.tight_layout()
plt.show()
Note that the default PSF_UPSCALE is 4.
Result of Adding an Artificial Source Tutorial
Comparison of the astronomical scene before and after adding the PSF as described in Tutorial #4.
For additional questions not answered in this article, please contact the Roman Help Desk at STScI.
References
In addition to this documentation, STIPS is described in the following references. Users are encouraged to cite this publication:
- STIPS Development Team et al 2024, "STIPS: The Nancy Grace Roman Space Telescope Imaging Product Simulator", PASP 136 124502