Estimating source photometry#

In addition to the theoretical predictions of locations, kinematics and stellar distributions of sources, you may also wish to estimate the photometry of sources in your simulation. This tutorial will show you how to do this!

Learning Goals#

By the end of this tutorial you should know how to:

  • Compute the photometry of sources in a Population

  • Apply dust map extinction to the photometry of sources

  • Use custom distances for non-Milky Way simulations

Beware - extra dependencies required here!

You’ll need to have installed the extra dependencies of cogsworth to estimate photometry of sources. Check out the installation page for more details on how to do this! (You’ll probably just need to run pip install cogsworth[extras])

[1]:

import cogsworth
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import astropy.units as u

[2]:

# this all just makes plots look nice
%config InlineBackend.figure_format = 'retina'

plt.rc('font', family='serif')
plt.rcParams['text.usetex'] = False
fs = 24

# update various fontsizes to match
params = {'figure.figsize': (12, 8),
          'legend.fontsize': fs,
          'axes.labelsize': fs,
          'xtick.labelsize': 0.9 * fs,
          'ytick.labelsize': 0.9 * fs,
          'axes.linewidth': 1.1,
          'xtick.major.size': 7,
          'xtick.minor.size': 4,
          'ytick.major.size': 7,
          'ytick.minor.size': 4}
plt.rcParams.update(params)
pd.options.display.max_columns = 999

Population photometry#

There are a couple of components that go into estimating the photometry of sources in a simulated Population:

  • Effective temperature, surface gravity and metallicity of each system (to get bolometric luminosity)

  • Distance to each source (to convert absolute to apparent magnitudes)

  • Dust extinction for each source (to apply to apparent magnitudes)

The information for the bolometric luminosity is contained entirely in each COSMIC output bpp table and is calculated using the MIST isochrones. For the dust extinction we use the Bayestar2019 dust map by default and the distances can either be inputted or calculated assuming a Milky Way setup and that the observer is on Earth.

Calculating photometry is fairly simple given an evolved population - let’s try this out.

[3]:

p = cogsworth.pop.Population(2500, use_default_BSE_settings=True)
p.create_population()

Run for 2500 binaries
Sampled 3348 binaries
[1e-02s] Sample initial binaries
[0.7s] Evolve binaries (run COSMIC)

Integrating orbits: 100%|██████████| 3356/3356 [00:04<00:00, 740.76it/s]

[6.0s] Integrate galactic orbits (run gala)
Overall: 6.8s

[4]:

p.get_observables(filters=["Gaia_G_EDR3", "Gaia_BP_EDR3", "Gaia_RP_EDR3"],
                  ignore_extinction=False, assume_mw_galactocentric=True,
                  silence_bounds_warning=True);

This populates the p.observables table with:

  • The assumed extinction for each source

  • Bolometric magnitudes

  • Magnitudes in each filter

  • For the G filter we also output whether the secondary star is brighter (and the observed temperature and surface gravity based on this)

[5]:

p.observables

[5]:
Av_1 Av_2 M_abs_1 m_app_1 M_abs_2 m_app_2 Gaia_G_EDR3_app_1 Gaia_G_EDR3_app_2 teff_obs log_g_obs secondary_brighter Gaia_G_EDR3_abs_1 Gaia_G_EDR3_abs_2 Gaia_BP_EDR3_app_1 Gaia_BP_EDR3_app_2 Gaia_BP_EDR3_abs_1 Gaia_BP_EDR3_abs_2 Gaia_RP_EDR3_app_1 Gaia_RP_EDR3_app_2 Gaia_RP_EDR3_abs_1 Gaia_RP_EDR3_abs_2
0 3.333000 6.0 10.273020 24.902600 11.526079 26.155659 27.786380 inf 3237.409615 4.972553 False 11.143872 inf 30.708783 inf 12.813383 inf 26.322171 inf 9.896075 inf
1 5.148000 6.0 9.205454 22.937900 9.818399 23.550844 26.270199 inf 3613.601960 4.930603 False 9.420474 inf 29.252073 inf 10.581822 inf 24.805402 inf 8.355534 inf
2 6.000000 6.0 8.145901 22.038654 9.927306 23.820058 25.977614 inf 3783.218240 4.783519 False 8.630722 inf 29.164091 inf 9.733847 inf 24.485698 inf 7.586910 inf
3 3.531000 6.0 9.034814 23.516264 9.175367 23.656817 26.155724 inf 3560.665698 4.889564 False 9.186757 inf 28.833339 inf 10.477703 inf 24.728263 inf 8.066396 inf
4 6.000000 6.0 8.854705 23.529008 10.092647 24.766949 27.398663 inf 3622.894135 4.876995 False 9.338264 inf 30.762635 inf 10.567622 inf 25.887019 inf 8.239965 inf
... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ...
3343 4.983000 6.0 7.376696 21.760561 11.144605 25.528469 25.230796 inf 4036.222484 4.671091 False 7.779006 inf 27.678958 inf 8.641159 inf 23.844498 inf 6.865615 inf
3344 3.892272 6.0 5.245240 19.613696 7.726504 22.094961 22.267580 inf 5440.884504 4.484475 False 5.094196 inf 23.747926 inf 5.495327 inf 21.088685 inf 4.520856 inf
3345 6.000000 6.0 6.553501 21.096736 10.502945 25.046180 25.025406 inf 4745.956610 4.644648 False 6.597684 inf 27.353759 inf 7.112288 inf 23.646858 inf 5.921349 inf
3346 6.000000 6.0 8.678926 23.146658 11.734802 26.202534 27.174966 inf 3742.484484 4.876034 False 9.177240 inf 30.142328 inf 10.161111 inf 25.706517 inf 8.191595 inf
3347 1.122000 6.0 9.119755 23.239831 12.498029 26.618105 24.660818 inf 3601.924721 4.914646 False 9.807927 inf 26.193570 inf 10.979074 inf 23.468550 inf 8.737116 inf

3348 rows × 21 columns

We can use cogsworth to immediately transform this into a colour-magnitude diagram coloured by stellar type to see how things look. You can see that there’s a clear main sequence with some evolved massive stars and lots of white dwarfs.

[6]:

cogsworth.plot.plot_cmd(p, m_filter="Gaia_G_EDR3", c_filter_1="Gaia_BP_EDR3", c_filter_2="Gaia_RP_EDR3", s=5);
../../_images/tutorials_observables_photometry_10_0.png

Custom distances#

When computing photometry you can either assume a Milky Way galactocentric frame (as above) to automatically work out the distance to each system, or you can manually specify the distance to your sources.

Let’s try this for an example where we also turn off extinction entirely.

[7]:

p.get_observables(filters=["Gaia_G_EDR3", "Gaia_BP_EDR3", "Gaia_RP_EDR3"], ignore_extinction=True,
                  assume_mw_galactocentric=False, distances=np.ones(len(p.final_pos)) * u.kpc,
                  silence_bounds_warning=True)

[7]:
Av_1 Av_2 M_abs_1 m_app_1 M_abs_2 m_app_2 Gaia_G_EDR3_app_1 Gaia_G_EDR3_app_2 teff_obs log_g_obs secondary_brighter Gaia_G_EDR3_abs_1 Gaia_G_EDR3_abs_2 Gaia_BP_EDR3_app_1 Gaia_BP_EDR3_app_2 Gaia_BP_EDR3_abs_1 Gaia_BP_EDR3_abs_2 Gaia_RP_EDR3_app_1 Gaia_RP_EDR3_app_2 Gaia_RP_EDR3_abs_1 Gaia_RP_EDR3_abs_2
0 0.0 0.0 10.273020 20.273020 11.526079 21.526079 21.143872 inf 3237.409615 4.972553 False 11.143872 inf 22.813383 inf 12.813383 inf 19.896075 inf 9.896075 inf
1 0.0 0.0 9.205454 19.205454 9.818399 19.818399 19.420474 inf 3613.601960 4.930603 False 9.420474 inf 20.581822 inf 10.581822 inf 18.355534 inf 8.355534 inf
2 0.0 0.0 8.145901 18.145901 9.927306 19.927306 18.630722 inf 3783.218240 4.783519 False 8.630722 inf 19.733847 inf 9.733847 inf 17.586910 inf 7.586910 inf
3 0.0 0.0 9.034814 19.034814 9.175367 19.175367 19.186757 inf 3560.665698 4.889564 False 9.186757 inf 20.477703 inf 10.477703 inf 18.066396 inf 8.066396 inf
4 0.0 0.0 8.854705 18.854705 10.092647 20.092647 19.338264 inf 3622.894135 4.876995 False 9.338264 inf 20.567622 inf 10.567622 inf 18.239965 inf 8.239965 inf
... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ...
3343 0.0 0.0 7.376696 17.376696 11.144605 21.144605 17.779006 inf 4036.222484 4.671091 False 7.779006 inf 18.641159 inf 8.641159 inf 16.865615 inf 6.865615 inf
3344 0.0 0.0 5.245240 15.245240 7.726504 17.726504 15.094196 inf 5440.884504 4.484475 False 5.094196 inf 15.495327 inf 5.495327 inf 14.520856 inf 4.520856 inf
3345 0.0 0.0 6.553501 16.553501 10.502945 20.502945 16.597684 inf 4745.956610 4.644648 False 6.597684 inf 17.112288 inf 7.112288 inf 15.921349 inf 5.921349 inf
3346 0.0 0.0 8.678926 18.678926 11.734802 21.734802 19.177240 inf 3742.484484 4.876034 False 9.177240 inf 20.161111 inf 10.161111 inf 18.191595 inf 8.191595 inf
3347 0.0 0.0 9.119755 19.119755 12.498029 22.498029 19.807927 inf 3601.924721 4.914646 False 9.807927 inf 20.979074 inf 10.979074 inf 18.737116 inf 8.737116 inf

3348 rows × 21 columns

You’ll notice this runs a lot faster now that we don’t need to query a dustmap or transform any coordinates!

Wrap-up#

Now you know all about how to transform your intrinsic population into an observable sample, nice job! Keep reading the next tutorial to learn about how you can use these observable to tell whether a given system is observable by Gaia.

Note

This tutorial was generated from a Jupyter notebook that can be found here.