๐ Quickstartยถ
This page shows how to load a time-series spectroscopic dataset, do some basic calculations with it, generate some visualizations, and then save it out to another format.
from chromatic import download_from_mast, read_rainbow, version
from astroquery.mast import Observations
import astropy.units as u
version()
'0.4.14'
๐พ Downloadยถ
Let's download the JWST Early Release Observation of HAT-P-18b, one of first exoplanet transit datasets to be gathered by the telescope. We'll get the default pipeline x1dints (Stage 2) outputs; there are lots of reasons why we shouldn't use these particular pipeline files for science, but they're useful for a quick initial look.
downloaded = download_from_mast(
proposal_id="2734", instrument_name="NIRISS/SOSS", target_name="WASP-96"
)
INFO: Found cached file ./mastDownload/JWST/jw02734002001_04101_00001-seg001_nis/jw02734002001_04101_00001-seg001_nis_x1dints.fits with expected size 55373760. [astroquery.query]
INFO: Found cached file ./mastDownload/JWST/jw02734002001_04101_00001-seg002_nis/jw02734002001_04101_00001-seg002_nis_x1dints.fits with expected size 55373760. [astroquery.query]
INFO: Found cached file ./mastDownload/JWST/jw02734002001_04101_00001-seg003_nis/jw02734002001_04101_00001-seg003_nis_x1dints.fits with expected size 44308800. [astroquery.query]
downloaded
| Local Path | Status | Message | URL |
|---|---|---|---|
| str106 | str8 | object | object |
| ./mastDownload/JWST/jw02734002001_04101_00001-seg001_nis/jw02734002001_04101_00001-seg001_nis_x1dints.fits | COMPLETE | None | None |
| ./mastDownload/JWST/jw02734002001_04101_00001-seg002_nis/jw02734002001_04101_00001-seg002_nis_x1dints.fits | COMPLETE | None | None |
| ./mastDownload/JWST/jw02734002001_04101_00001-seg003_nis/jw02734002001_04101_00001-seg003_nis_x1dints.fits | COMPLETE | None | None |
๐งโ๐ป Readยถ
Next, let's load that transit dataset into a Rainbow (๐) object. These chromatic ๐ objects keep track of how the brightness of source changes across both wavelength and time.
filenames = downloaded["Local Path"]
rainbow = read_rainbow(filenames)
๐๐ค This file contains data for 2 spectroscopic orders. Because no `order=` keyword was supplied, we're defaulting to first order. You can hide this warning by expliciting stating which order you want to load. For this file, the options include [1 2].
๐๐ค The 2048 input wavelengths were not monotonically increasing. <๐(2048w, 280t)> has been sorted from lowest to highest wavelength. If you want to recover the original wavelength order, the original wavelength indices are available in `rainbow.original_wave_index`.
The ๐ object we just loaded provides easy access to the different dimensions we might want from the dataset, arrays like wavelength, time, flux, or uncertainty. If appropriate, quantities will have astropy Units.
rainbow.wavelength
rainbow.time
rainbow.flux
rainbow.uncertainty
๐งฎ Calculateยถ
The absolute flux doesn't matter that much for many transit analyses, so let's use .normalize() to normalize out the median spectrum of the star, converting the data to relative brightness within each wavelength. This ๐ action returns another ๐ object, just with the brightness normalized.
normalized = rainbow.normalize()
The dataset is really large, so for making some simple visualization it might help to average over bins of wavelength and/or time. Let's use .bin() to bin onto a (logarithmically) uniform wavelength grid, returning the binned ๐.
binned = normalized.bin(R=200, dt=4 * u.minute)
/Users/zach/opt/anaconda3/envs/2024/lib/python3.12/site-packages/astropy/units/quantity.py:1865: RuntimeWarning: All-NaN slice encountered result = super().__array_function__(function, types, args, kwargs)
The times in this dataset are measured relative to some arbitrary time in the distant past. To make them easier to interpret we can phase-fold the times so they're measured relative to the mid-transit time, when the planet is directly between the star and us, according to the planet's orbital properties from the NASA Exoplanet Archive. Let's use .fold() to change the times, returning a phase-folded ๐.
folded = binned.fold(period=5.5080232 * u.day, t0=2459033.31735 * u.day)
๐จ Visualizeยถ
We can visualize the dataset by making a map of the star's brightness across both wavelength and time, an image in which each the brightness along row corresponds to a transit light curve at that wavelength. Let's use .imshow() to create this map for the normalized, binned, folded ๐.
folded.imshow();
It might be nice to look closely at the light curves within a particular wavelength range. Let's use .imshow_interact() to interactively explore the ๐. Click and drag on the panel on the left to select the wavelegnth range to display as a light curve on the right.
folded.imshow_interact()
What a dataset! It looks like there's something a little odd happening at about 2 microns (probably contamination from another star or spectrograph order) and a starspot crossing just after mid-transit, but otherwise it's a remarkably beautiful transit from a very impressive telescope!
๐งถ Buildยถ
Because many of the actions possible with ๐ objects return other ๐ objects, it's possible to connect multiple steps into a single command, building up complicated analysis stories with relatively succinct code.
(
rainbow.normalize()
.flag_outliers()
.bin(R=10)
.fold(period=5.5080232 * u.day, t0=2459033.31735 * u.day)
.plot()
);
/Users/zach/opt/anaconda3/envs/2024/lib/python3.12/site-packages/astropy/units/quantity.py:1865: RuntimeWarning: All-NaN slice encountered result = super().__array_function__(function, types, args, kwargs)
๐พ Saveยถ
Let's convert these data into a different format by saving it as a new file, which we might send around to share with our colleagues or publish along with a paper. chromatic can read and save ๐ datasets with a variety of formats, to try to ease collaboration across different pipelines and toolkits.
rainbow.save("jwst-hatp18b.rainbow.npy")
๐ Learnยถ
That's it! This quick tutorial highlighted chromatic's abilities to...
- load in time-series spectra or multiwavelength light curves from formats like
x1dints - access core data variables like
wavelength,time,flux,uncertainty - perform calculations like
.normalize,.bin,.fold - visualize the data with
.imshow,.imshow_interact,.plot
Hopefully, you're now curious to read through the User Guide to learn more about options for reading ๐s, doing actions with ๐s, visualizing ๐s in different ways, and more! You can also run the .help() method associated with any ๐ object to get a quick summary of what other methods are available for it>
rainbow.help()
Hooray for you! You asked for help on what you can do with this ๐ object. Here's a quick reference of a few available options for things to try. ----------- | actions | ----------- ๐๐งฎ๐ | +-*/ Do basic math operations with two Rainbows. ๐๐๐ช | .[:,:]() Index, slice, or mask a Rainbow to get a subset. ๐๐ง๐ | .align_wavelengths() Align spectra with different wavelength grids onto one shared axis. ๐๐งบ๐งฑ | .bin() Bin to a new wavelength or time grid. ๐๐งโ๐คโ๐ง๐ | .compare() Connect to other ๐ objects for easy comparison. ๐๐โฐ | .concatenate_in_time() Stitch together two Rainbows with identical wavelengths. ๐๐๐ | .concatenate_in_wavelength() Stitch together two Rainbows with identical times. ๐๐ฉ๐ | .flag_outliers() Flag outlier data points. ๐โฒ๐ | .fold() Fold times relative to a particular period and epoch. ๐๐งบโฐ | .get_average_lightcurve_as_rainbow() Bin down to a single integrated light curve. ๐๐งบ๐ | .get_average_spectrum_as_rainbow() Bin down to a single integrated spectrum. ๐๐ง๐ฒ | .inject_noise() Inject (uncorrelated, simple) random noise. ๐๐ง๐น | .inject_systematics() Inject (correlated, wobbly) systematic noise. ๐โญ๏ธ๐ป | .inject_spectrum() Inject a static stellar spectrum. ๐๐ช๐ | .inject_transit() Inject a transit signal. ๐๐ซ๐ | .normalize() Normalize by dividing through by a typical spectrum (and/or light curve). ๐๐ดโโ ๏ธ๐ | .remove_trends() Remove smooth trends in time and/or wavelength. ๐๐๐ | .shift() Doppler shift wavelengths. ๐๐ฑ๐ | .trim() Trim away wavelengths or times. ---------------- | get/timelike | ---------------- โฐ๐๐ | .get_average_lightcurve() Get the weighted average light curve. โฐ๐๐ | .get_for_time() Get a quantity associated with a time index. โฐ๐๐ | .get_median_lightcurve() Get the median light curve. โฐ๐๐ | .get_ok_data_for_time() Get a quantity associated with a time index. โฐ๐ฐ๐ | .get_times_as_astropy() Get the times as an astropy Time object. โฐ๐๐ | .set_times_from_astropy() Set the times from an astropy Time object (modifies in-place). ---------------- | get/wavelike | ---------------- ๐๐๐ | .get_average_spectrum() Get the weighted average spectrum. ๐๐๐ | .get_for_wavelength() Get a quantity associated with a wavelength index. ๐๐ฏ๐ | .get_measured_scatter() Get the measured scatter on the time series for each wavelength. ๐๐๐ | .get_median_spectrum() Get the median spectrum. ๐๐๐ | .get_ok_data_for_wavelength() Get a quantity associated with a wavelength index. ๐๐๐ | .get_spectral_resolution() Get the spectral resolution (R=w/dw). ----------- | helpers | ----------- ๐๐๐ | .help() Get one-line help summaries of available methods. ๐๐๐ชต | .history() Get the history that went into this Rainbow. ๐พ๐๐ผ | .save() Save this Rainbow out to a permanent file. ------------------ | visualizations | ------------------ ๐จ๐ฝโฐ | .animate_lightcurves() Animate a sequence of light curves across different wavelengths. ๐จ๐ฝ๐ | .animate_spectra() Animate a sequence of spectra across different times. ๐จ๐ผ๐บ | .imshow() Paint a map of flux across wavelength and time. ๐จ๐น๐บ | .imshow_interact() Show flux map and lightcurves with interactive wavelength selection. ๐จ๐๐บ | .pcolormesh() Paint a map of flux across wavelength and time (with non-uniform grids). ๐จ๐๐งถ | .plot() Plot a sequence of light curves with vertical offsets. ------------------------------ | visualizations/diagnostics | ------------------------------ ๐จ๐๐บ | .imshow_quantities() Show multiple 2D (wavelength and time) quantities as imshow maps. ๐จ๐๐งถ | .plot_quantities() Show multiple 1D (wavelength or time) quantities as scatter plots. --------------------------- | visualizations/timelike | --------------------------- ๐จโฐ๐ | .plot_average_lightcurve() Plot the weighted average flux per time. ๐จโฐ๐ | .plot_median_lightcurve() Plot the median flux per time. --------------------------- | visualizations/wavelike | --------------------------- ๐จ๐๐ญ | .plot_average_spectrum() Plot the weighted average flux per wavelength. ๐จ๐๐ | .plot_spectral_resolution() Plot the spectral resolution per wavelength. ๐จ๐๐ง | .plot_noise_comparison() Plot the measured and expected scatter per wavelength. ๐จ๐๐ฅฆ | .plot_noise_comparison_in_bins() Plot measured and expected scatter in different size bins.