Simulating from Multiple Surveys
LightCurveLynx can simulate observations from multiple surveys in a single run. In this notebook, we show how to perform these simulations and discuss how it works.
[1]:
import numpy as np
from lightcurvelynx.astro_utils.passbands import PassbandGroup
from lightcurvelynx.models.basic_models import ConstantSEDModel
from lightcurvelynx.obstable.opsim import OpSim
from lightcurvelynx.simulate import simulate_lightcurves
from lightcurvelynx.survey_info import SurveyInfo
# Usually we would not hardcode the path to the passband files, but for this demo we will use a relative path
# to the test data directory so that we do not have to download the files.
table_dir = "../../tests/lightcurvelynx/data/passbands"
/home/docs/checkouts/readthedocs.org/user_builds/lightcurvelynx/envs/latest/lib/python3.12/site-packages/tqdm/auto.py:21: TqdmWarning: IProgress not found. Please update jupyter and ipywidgets. See https://ipywidgets.readthedocs.io/en/stable/user_install.html
from .autonotebook import tqdm as notebook_tqdm
Prerequisite Data
For each survey that we would like to include in the simulation, we need two pieces of information:
An
ObsTablethat includes the survey’s pointing and noise information.A
PassbandGroupfor that survey.
Let’s start by creating toy data for each survey.
The first survey will include pointings at two locations (0.0, 10.0) and (180.0, -10.0) in the “g” and “r” bands.
[2]:
obsdata1 = {
"time": [0.0, 1.0, 2.0, 3.0],
"ra": [0.0, 0.0, 180.0, 180.0],
"dec": [10.0, 10.0, -10.0, -10.0],
"filter": ["g", "r", "g", "r"],
"zp": [5.0, 6.0, 7.0, 8.0],
"seeing": [1.12, 1.12, 1.12, 1.12],
"skybrightness": [20.0, 20.0, 20.0, 20.0],
"exptime": [29.2, 29.2, 29.2, 29.2],
"nexposure": [2, 2, 2, 2],
}
obstable1 = OpSim(obsdata1)
passband_group1 = PassbandGroup.from_preset(
preset="LSST",
table_dir=table_dir,
filters=["g", "r", "i"],
)
survey_info1 = SurveyInfo(obstable=obstable1, passbands=passband_group1, survey_name="LSST")
The second survey will include pointings at two locations (0.0, 10.0) and (90.0, -10.0) in the “r” and “z” bands.
[3]:
obsdata2 = {
"time": [0.5, 1.5, 2.5, 3.5],
"ra": [0.0, 90.0, 0.0, 90.0],
"dec": [10.0, -10.0, 10.0, -10.0],
"filter": ["r", "i", "r", "i"],
"zp": [1.0, 2.0, 3.0, 4.0],
"seeing": [1.12, 1.12, 1.12, 1.12],
"skybrightness": [20.0, 20.0, 20.0, 20.0],
"exptime": [29.2, 29.2, 29.2, 29.2],
"nexposure": [2, 2, 2, 2],
}
obstable2 = OpSim(obsdata2)
passband_group2 = PassbandGroup.from_preset(
preset="LSST",
table_dir=table_dir,
filters=["r", "i"],
)
survey_info2 = SurveyInfo(obstable=obstable2, passbands=passband_group2, survey_name="LSST")
Note that, while we use the OpSim objects for both ObsTables and we use LSST preset for both passbands, this is only for simplicity of the example (only some of the data files are installed by default).
In most cases users will be interested in combining data from multiple surveys. For example we could simulate a trio of surveys: Rubin ( OpSim and LSST passbands), ZTF (ZTFObsTable with ZTF passbands) and Roman (data coming soon).
Model Creation
Nothing changes in the model creation step. We define a model and its parameters as we would with any other simulation. The model is just a recipe for generating the noise-free fluxes given parameter values. It doesn’t need to know anything about the survey.
Here we use a constant SED model (same value for all times and wavelengths). We place the object at (0.0, 10.0) so it is observed by some of the pointings from each survey.
[4]:
model = ConstantSEDModel(brightness=100.0, t0=0.0, ra=0.0, dec=10.0, redshift=0.0, node_label="my_star")
If we want to sample the positions (RA, dec) from the footprints of the combined surveys, we can use the ApproximateMOCSampler. That class contains built-in constructor methods for creating a sampling footprint from either the union or intersection of footprints of the individual surveys. See the sampling_positions notebook for details.
Simulation
We run the simulation by passing in a list of SurveyInfo (note this changed in v0.4.0. In previous versions, users passed lists of ObsTable and PassbandGroup). The parameter space is only sampled once for each object, so the observations in each survey are consistent with respect to the parameterization. The times of observation and filters used are determined by each survey. And the bandflux is computed using that survey’s passbands.
We can simulate a single object and look at its light curve.
[5]:
results = simulate_lightcurves(
model,
1,
[survey_info1, survey_info2],
obstable_save_cols=["zp"],
)
print(results["lightcurve"][0])
Simulating: 100%|██████████| 1/1 [00:00<00:00, 310.53obj/s]
mjd filter flux fluxerr flux_perfect survey_idx obs_idx \
0 0.0 g -728.199199 895.605471 100.0 0 0
1 0.5 r -342.666791 338.724102 100.0 1 0
2 1.0 r 428.160170 1015.412583 100.0 0 1
3 2.5 r 173.489453 642.443596 100.0 1 2
is_saturated zp
0 False 5.0
1 False 1.0
2 False 6.0
3 False 3.0
In addition to the standard information like mjd and flux, we capture the survey index to help us trace which survey the data came from. This is stored in the survey_idx column of the nested light curve table. We also can save per-survey, per-observation data, such as the zero point for that observation, using the obstable_save_cols parameter as shown above.
Adding Spectra
The list of “PassbandGroups” can also include Spectrograph objects in order to evalue spectra. Let’s add a third survey with a spectrograph. Since spectographic surveys use different information, we can use the simpler SpectrographObsTable.
[6]:
from lightcurvelynx.astro_utils.spectrograph import Spectrograph
from lightcurvelynx.obstable.spectrograph_table import SpectrographObsTable
obsdata3 = {
"time": [0.0, 1.0, 2.0, 3.0],
"ra": [0.0, 0.0, 180.0, 180.0],
"dec": [10.0, 10.0, -10.0, -10.0],
}
obstable3 = SpectrographObsTable(obsdata3)
spectograph = Spectrograph(np.arange(4000, 8000, 100))
survey_info3 = SurveyInfo(obstable=obstable3, passbands=spectograph, survey_name="spectrograph")
We add the spectrograph information as another survey and “passbandgroup”.
[7]:
results2 = simulate_lightcurves(
model,
1,
[survey_info1, survey_info2, survey_info3],
obstable_save_cols=["zp"],
)
results2.head()
Simulating: 100%|██████████| 1/1 [00:00<00:00, 337.11obj/s]
[7]:
| id | ra | dec | nobs | t0 | z | lightcurve | params | spectra | ||||||||||||||||||||||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 0 | 0.000000 | 10.000000 | 6 | 0.000000 | 0.000000 |
|
{'my_star.ra': 0.0, 'my_star.dec': 10.0, 'my_star.redshift': 0.0, 'my_star.t0': 0.0, 'my_star.distance': None, 'my_star.brightness': 100.0} |
|
The results now have two nested columns. As before, the ‘lightcurve’ column contains all of the light curve points from all of the (non-spectrographic) surveys. A new ‘spectrograph’ column contains each recorded spectra.
FAQ
Can we use multiple surveys and parallelization?
Yes.
What happens if the filter names overlap?
For the simulation, everything will just work. Since you are passing one PassbandGroup for each ObsTable it will only look for the matching filter names in that group. So if you are using LSST data with observations in r and then ZTF data with observations in r, the SEDs will be integrated by the correct filter transmission at each step.
However the name of the filter in the result column will match what appears in the ObsTable columns. So you may only see ‘r’ for each entry. If you are interested in plotting the light curves or performing any other per-filter analysis, you can expand the filter names with the survey prefix using test_results_use_full_filter_names()
[8]:
from lightcurvelynx.utils.post_process_results import results_use_full_filter_names
print("Filters (before):", list(results["lightcurve.filter"].values))
results_use_full_filter_names(results, [passband_group1, passband_group2])
print("Filters (after):", list(results["lightcurve.filter"].values))
Filters (before): ['g', 'r', 'r', 'r']
Filters (after): ['LSST_g', 'LSST_r', 'LSST_r', 'LSST_r']