WFC3/UVIS Filter Transformations with stsynphot

---

Learning Goals

By the end of this tutorial, you will:

• Generate synthetic observations using synphot and stsynphot.

• Find color terms between WFC3/UVIS filters and non-HST filters.

• Plot bandpasses to investigate various throughputs.

Table of Contents

Introduction
1. Imports
2. Select filters for the transformation
3. Define a spectrum
4. Select UVIS chips
5. Select magnitude systems
6. Generate outputs
7. Plot bandpasses
8. Conclusions
Additional Resources
About the Notebook
Citations

Introduction

This notebook computes color terms between selected WFC3/UVIS filters and non-HST filters, such as Johnson-Cousins, for a user-defined reference spectrum. The terms as given are the difference between the magnitude of the spectrum in the selected non-HST filter and the corresponding UVIS filter.

This tool reproduces the methods described in section 4 of WFC3 ISR 2014-16, but will automatically use the latest available spectra and throughput tables.

stsynphot requires access to data distributed by the Calibration Data Reference System (CRDS) in order to operate. Both packages look for an environment variable called PYSYN_CDBS to find the directory containing these data.

Users can obtain these data files from the CDRS. Information on how to obtain the most up-to-date reference files (and what they contain) can be found here. An example of how to download the files with curl and set up this environment variable is presented below.

For detailed instructions on how to install and set up these packages, see the synphot and stsynphot documentation.

1. Imports

This notebook assumes you have created the virtual environment in WFC3 notebooks’ installation instructions.

We import:

• os for setting environment variables

• tarfile for extracting a .tar archive

• numpy for handling array functions

• pandas for managing data

• matplotlib.pyplot for plotting data

• astropy.units and synphot.units for handling units

• synphot and stsynphot for evaluating synthetic photometry

Additionally, we will need to set the PYSYN_CDBS environment variable before importing stsynphot. We will also create a Vega spectrum using synphot’s inbuilt from_vega() method, as the latter package will supercede this method’s functionality and require a downloaded copy of the latest Vega spectrum to be provided.

import os
import tarfile
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import astropy.units as u
import synphot.units as su
import synphot as syn
from synphot import Observation

%matplotlib inline

vegaspec = syn.SourceSpectrum.from_vega()

This section obtains the WFC3 throughput component tables for use with stsynphot. This step only needs to be done once. If these reference files have already been downloaded, this section can be skipped.

!curl -O https://archive.stsci.edu/hlsps/reference-atlases/hlsp_reference-atlases_hst_multi_everything_multi_v11_sed.tar

  % Total    % Received % Xferd  Average Speed   Time    Time     Time  Current
                                 Dload  Upload   Total   Spent    Left  Speed

  0     0    0     0    0     0      0      0 --:--:-- --:--:-- --:--:--     0
  0     0    0     0    0     0      0      0 --:--:-- --:--:-- --:--:--     0

  0  796M    0 8112k    0     0  7906k      0  0:01:43  0:00:01  0:01:42 7898k

  6  796M    6 48.1M    0     0  23.9M      0  0:00:33  0:00:02  0:00:31 23.9M

 11  796M   11 90.3M    0     0  30.0M      0  0:00:26  0:00:03  0:00:23 30.0M

 16  796M   16  132M    0     0  33.0M      0  0:00:24  0:00:04  0:00:20 33.0M

 21  796M   21  173M    0     0  34.5M      0  0:00:23  0:00:05  0:00:18 34.6M

 26  796M   26  214M    0     0  35.5M      0  0:00:22  0:00:06  0:00:16 41.2M

 32  796M   32  255M    0     0  36.4M      0  0:00:21  0:00:07  0:00:14 41.4M

 37  796M   37  296M    0     0  36.9M      0  0:00:21  0:00:08  0:00:13 41.1M

 42  796M   42  337M    0     0  37.5M      0  0:00:21  0:00:09  0:00:12 41.0M

 47  796M   47  378M    0     0  37.8M      0  0:00:21  0:00:10  0:00:11 41.0M

 52  796M   52  420M    0     0  38.0M      0  0:00:20  0:00:11  0:00:09 41.0M

 57  796M   57  460M    0     0  38.3M      0  0:00:20  0:00:12  0:00:08 40.9M

 63  796M   63  501M    0     0  38.5M      0  0:00:20  0:00:13  0:00:07 41.0M

 68  796M   68  543M    0     0  38.7M      0  0:00:20  0:00:14  0:00:06 41.0M

 73  796M   73  584M    0     0  38.9M      0  0:00:20  0:00:15  0:00:05 41.1M

 78  796M   78  626M    0     0  39.0M      0  0:00:20  0:00:16  0:00:04 41.3M

 83  796M   83  667M    0     0  39.2M      0  0:00:20  0:00:17  0:00:03 41.4M

 88  796M   88  708M    0     0  39.2M      0  0:00:20  0:00:18  0:00:02 41.1M

 94  796M   94  748M    0     0  39.3M      0  0:00:20  0:00:19  0:00:01 41.0M

 99  796M   99  790M    0     0  39.5M      0  0:00:20  0:00:20 --:--:-- 41.2M

100  796M  100  796M    0     0  39.5M      0  0:00:20  0:00:20 --:--:-- 41.2M

Once the downloaded is complete, extract the file and set the environment variable PYSYN_CDBS to the path of the trds subdirectory. The next cell will do this for you, as long as the .tar file downloaded above has not been moved.

tar_archive = 'hlsp_reference-atlases_hst_multi_everything_multi_v11_sed.tar'
extract_to = 'hlsp_reference-atlases_hst_multi_everything_multi_v11_sed'
abs_extract_to = os.path.abspath(extract_to)

with tarfile.open(tar_archive, 'r') as tar:
    for member in tar.getmembers():
        member_path = os.path.abspath(os.path.join(abs_extract_to, member.name))
        if member_path.startswith(abs_extract_to):
            tar.extract(member, path=extract_to)   
        else:
            print(f"Skipped {member.name} due to potential security risk")

os.environ['PYSYN_CDBS'] = os.path.join(abs_extract_to, 'grp/redcat/trds/')

Now, after having set up PYSYN_CDBS, we import stsynphot. A warning regarding the Vega spectrum is expected here.

import stsynphot as stsyn

WARNING: Failed to load Vega spectrum from /home/runner/work/hst_notebooks/hst_notebooks/notebooks/WFC3/filter_transformations/hlsp_reference-atlases_hst_multi_everything_multi_v11_sed/grp/redcat/trds//calspec/alpha_lyr_stis_011.fits; Functionality involving Vega will be severely limited: FileNotFoundError(2, 'No such file or directory') [stsynphot.spectrum]

2. Select filters for the transformation

Define the filters to use for computing the transformation. One filter should be a UVIS filter, and the other a non-HST filter such as a Johnson-Cousins filter.

Filter names should be input as a list of tupled strings. Each tuple represents a pair of filters to convert between, and should contain the non-HST filter as the first element, and the UVIS filter as the second.

For non-HST filters, the filter system be included in the string, separated from the filter name by a comma (e.g. 'johnson, v' or 'sdss, g'). The available non-HST filters are listed here:

System	Bands
cousins	r, i
galex	nuv, fuv
johnson	u, b, v, r, i, j, k
landolt	u, b, v, r, i
sdss	u, g, r, i, z,
stromgren	u, v, b, y

Furthermore, Johnson-Cousins filters with corresponding UVIS filters are listed here:

Johnson-Cousins Filter	UVIS Filter
U	F336W
B	F475W
V	F555W/F606W
I	F814W

A summary of the UVIS filters, with descriptions, is available in Section 6.5.1 of the WFC3 Instrument Handbook

The notebook is currently set up to return the color terms between the V and I Johnson-Cousins filters, and corresponding UVIS filters.

filter_pairs = [('johnson, v', 'f555w'), ('cousins, i', 'f814w')]

3. Define a spectrum

Define a spectrum to get color terms for. Some common options are embedded below. A wide array of reference spectra are available for download from spectral atlases located here.

# Blackbody (5000 K)
blackbody_temperature = 5000

model = syn.models.BlackBody1D(blackbody_temperature)
source_spectrum  = syn.SourceSpectrum(model)

# Power law 
pl_index = 0

model = syn.models.PowerLawFlux1D(amplitude=flux_in, x_0=wl_in, alpha=pl_index)
source_spectrum  = syn.SourceSpectrum(model)
                   
# Load from a FITS table (e.g. a CALSPEC spectrum)
source_spectrum = syn.SourceSpectrum.from_file('/path/to/your/spectrum.fits')

Currently, the notebook is configured to use a 5000 K blackbody spectrum.

blackbody_temperature = 5000

model = syn.models.BlackBody1D(blackbody_temperature)
source_spectrum = syn.SourceSpectrum(model)

4. Select UVIS chips

Quantum efficiency differences between the two UVIS chips mean that you must specify which chips to use for computing color terms. Simply set the chip you would like to use to True and the other to False, or set both to True if you would like coefficients for both.

use_uvis1 = True
use_uvis2 = True
chips = [use_uvis1, use_uvis2]

5. Select magnitude systems

Select which magnitude systems you would like color terms to be provided for. Set those you would like to True and others to False.

ABMAG = True
STMAG = True
VEGAMAG = False

mags = [('ABMAG', u.ABmag, ABMAG), ('STMAG', u.STmag, STMAG),
        ('VEGAMAG', su.VEGAMAG, VEGAMAG)]

6. Generate outputs

Generate a data frame containing the color terms for the inputs you have specified.

First, let’s define the column names for the output table, and a list to fill with table rows.

cols = ['Filter Pair', 'Chip']
rows = []

Then, append the names of magnitude systems being used.

for name, _, toggle in mags:
    if toggle:
        cols.append(f'{name} Color Term')

Next, iterate over filter pairs. For each filter pair, this loop will:

• generate observation mode strings, bandpasses, and observations

• calculate the color term and append it

• append filters, chip, and color term as a row to rows

for pair in filter_pairs:
    comparison_filter, uvis_filter = pair  # Unpack filters
    filt_str = comparison_filter + ' - ' + uvis_filter

    for i, toggle in enumerate(chips):
        if not toggle:
            continue
        chip_str = 'UVIS' + str(i + 1)

        # Generate observation mode strings, bandpasses, observations
        comparison_obsmode = comparison_filter
        uvis_obsmode = 'wfc3, ' + chip_str + ', ' + uvis_filter

        comparison_bp = stsyn.band(comparison_obsmode)
        uvis_bp = stsyn.band(uvis_obsmode)

        comparison_observation = Observation(source_spectrum, comparison_bp)
        uvis_observation = Observation(source_spectrum, uvis_bp)
        row = [filt_str, chip_str]  # Append filters and chip to row

        for name, unit, toggle in mags:
            if not toggle:
                continue
            comparison_countrate = comparison_observation.effstim(
                flux_unit=unit, vegaspec=vegaspec)
            uvis_countrate = uvis_observation.effstim(
                flux_unit=unit, vegaspec=vegaspec)
            color = comparison_countrate - uvis_countrate  # Find color term
            row.append(f'{color.value:.3f}')  # Append color term

        rows.append(row)  # Append row to list of rows

Finally, generate and return the output table.

df = pd.DataFrame(rows, columns=cols)
df

	Filter Pair	Chip	ABMAG Color Term	STMAG Color Term
0	johnson, v - f555w	UVIS1	-0.094	-0.025
1	johnson, v - f555w	UVIS2	-0.094	-0.025
2	cousins, i - f814w	UVIS1	0.005	-0.039
3	cousins, i - f814w	UVIS2	0.005	-0.038

If you wish to save the output table as a .txt file, please uncomment and execute the code block below.

# df.to_csv('your/path/here.txt', sep='\t')

7. Plot bandpasses

It can be nice to see your selected bandpass pairs plotted with each other. The cell below will generate a figure with subplots for each filter pair specified above, as well as the relevant portion of the spectrum you’ve defined, all normalized to fit on the same axes.

Note: For the purposes of these plots, the non-HST bandpass and spectrum have been scaled to the amplitude of the HST bandpass, which reflects the actual total system throughput as a function of wavelength.

fig, axs = plt.subplots(1, len(filter_pairs), sharey=True, figsize=(
    4*len(filter_pairs), 5))  # Instantiate subplots
axs[0].set_ylabel('Throughput')
for i, pair in enumerate(filter_pairs):

    f1, f2 = pair
    bp1 = stsyn.band(f1)
    bp2 = stsyn.band('wfc3, uvis1,' + f2)

    # Create wavelength array for subplot based on average bandpass wavelength and width
    avgwave = (bp1.avgwave().to(u.nm) + bp2.avgwave().to(u.nm))/2
    width = (bp1.rectwidth().to(u.nm) + bp2.rectwidth().to(u.nm))/2
    left = max((avgwave - 1.5 * width).value, 1)
    right = (avgwave + 1.5 * width).value

    wl = np.arange(left, right) * u.nm

    # Normalize curves to fit on one set of axes
    bp1_norm = bp1(wl) / np.max(bp1(wl)) * np.max(bp2(wl))
    spec_norm = source_spectrum(
        wl) / np.max(source_spectrum(wl)) * np.max(bp2(wl))

    # Plot bandpasses and spectrum on subplot
    axs[i].plot(wl, bp1_norm, ls='--', label=f1, c='tab:blue')
    axs[i].plot(wl, bp2(wl),  ls='-.', label=f2, c='tab:red')
    axs[i].plot(wl, spec_norm, label='source spectrum', c='tab:purple')
    axs[i].set_xlabel('Wavelength (nm)')
    axs[i].legend(fontsize='x-small', loc='upper right')
plt.tight_layout()

8. Conclusions

Thank you for walking through this notebook. Now using WFC3 data, you should be more familiar with:

• Generating synthetic observations using synphot and stsynphot.

• Finding color terms between WFC3/UVIS filters and non-HST filters.

• Ploting bandpasses to investigate various throughputs.

Congratulations, you have completed the notebook!

Additional Resources

Below are some additional resources that may be helpful. Please send any questions through the HST Helpdesk.

• WFC3 Website

• WFC3 Instrument Handbook

• WFC3 Data Handbook

  • see sections 9.5.2 for reference to this notebook

About this Notebook

Authors: Aidan Pidgeon, Jennifer Mack; WFC3 Instrument Team

Updated on: 2021-09-13

Citations

If you use numpy, astropy, synphot, or stsynphot for published research, please cite the authors. Follow these links for more information about citing the libraries below:

• Citing numpy

• Citing astropy

• Citing synphot

• Citing stsynphot

---

Top of Page