Week 12: astroplanning

What I completed this week:

  • Submitted my synphot tutorials as a pull request to astropy/astropy-tutorials. Just waiting for someone (outside of my GSoC mentors) to review!
  • Moved my synphot work (including the wrapper and signal-to-noise functions) into an experimental folder. I’m still uncertain whether this folder will be merged as an experimental folder or if its contents will be merged back into the synphot API, but I guess we will see..! Also, I’m still open to suggestions on any other functionality I should include.
  • Emailed ESO inquiring about any plans to merge skycalc_cli into astroquery, and if they would be open to me submitting a PR (no response yet).
  • Thought more about how my wrapper, signal-to-noise functions, and “synspec” can be implemented into already-existing astronomical observing tools. Simply put, I think my GSoC work (and synphot in general) would serve as valuable additions to astroplan. As it stands, astroplan is a great tool for scheduling observations, but its definition leaves its functionality open to aspects of planning which are outside of the scheduling realm:

The goal of astroplan is to make a flexible toolbox for observation planning and scheduling. 

The astroplan Documentation

Specifically, I think a count rate function and exposure time calculator would be of greatest use. Then not only could astroplan help determine when to observe, it could also suggest exposure durations given model photometric/spectroscopic observations generated by synphot/synphot tools.

I also took a weekend trip with my family to go rafting on the Yakima river (pictured here). And I fed a cow. Its tongue was disturbingly pokey.

This week’s goals:

  • If everyone is open to the idea, I would like to submit a PR to astroplan for adding the exposure time calculator (supported by synphot and its tools)
Advertisements

Week 11: A wrapper says its first word

What I completed this week:

  • Continued working on the synphot wrapper class, which we’ve renamed “Exposure” (since “Observe” is a little too close to the other synphot class “Observation”). So far this class works pretty much like how I imagined it last week, where the most simple (read: default) synthetic observation can be initiated given a source spectrum, telescope diameter, and exposure time:
>>> from synphot.spectrum import SourceSpectrum
>>> from synphot.exposure import Exposure
>>> import astropy.units as u

>>> source_spec = SourceSpectrum(...)
>>> tele_diameter = 1 * u.m
>>> exptime = 1 * u.s

>>> some_observation = Exposure(source_spec, tele_diameter, exptime)

Then attributes can be called and methods applied to view some characteristics of this observation, such as the count rate, signal to noise ratio (SNR), and the exposure time needed to obtain a specified SNR:

>>> some_observation.aperture_area
<Quantity 0.78539816 m2>

>>> some_observation.countrate
<Quantity 126398.52488845 ct / s>

>>> some_observation.snr()
<Quantity 355.52570215>

>>> desired_snr = 100
>>> some_observation.required_exptime(desired_snr)
<Quantity 0.07911485 s>
I also took a little break to see family and go “tubing”, as seen above (wee)

This week’s goals:

  • Finish the wrapper! Now taking suggestions for anything else I should include
  • Take a look at telescopy and astroplan and think how they, along with synspec and any higher order functions we come up with, can be merged (okay, I really need to do this one this week)

Longer term goals:

Week 10: A wrapper is born

What I completed this week:

  • Continued working on the signal to noise functions pull request to synphot
    • Changed function names from howell_snr to ccd_snr and exptime_from_howell_snr to exptime_from_ccd_snr, which are more descriptive/informative
    • The units returned by ccd_snr are now dimensionless
  • Began writing a synphot wrapper. This work in progress PR is a bit of an embarrassing mess right now, but I’m envisioning something where you can initiate a class with the source spectrum of an object, give a telescope diameter, and boom, you can get attributes like SNR, countrate, exptime, etc… :
>>> from synphot.spectrum import SourceSpectrum
>>> from synphot.wrapper import Observe
>>> import astropy.units as u

>>> tele_diameter = 1 * u.m
>>> source_spec = SourceSpectrum(...)
>>> some_observation = Observe(source_spec, tele_diameter)

>>> some_observation.countrate
1e6 ct/s
>>> some_observation.snr
100
...

The observation can be customized by setting a bandpass, giving a desired SNR (to get a required exposure time), providing the quantum efficiency of the CCD or an atmospheric transmission model, etc.

But I’m not really sure if attributes are the best way to get these, er, “attributes” of the observation. Maybe argument-less methods would be better, e.g. some_observation.snr()? That might make it easier to alter some parameters like exposure time without having to create a new object… Anyway, this week has been full of thoughts like these.

How I kind of feel writing this. Something good will come together soon, I’m sure

This week’s goals:

  • Take a look at telescopy and astroplan and think how they, along with synspec and any higher order functions we come up with, can be merged

Longer term goals:

Week 9: synspec and other thoughts

What I completed this week:

  • Traveled to Seattle, where I’ll be for the rest of the summer
  • Added documentation on the signal to noise ratio (SNR) functions to the synphot docs
  • Integrated the SNR functions into the synphot tutorials I’ve been working on
  • Checked that the binned photometry (aka spectroscopy) is functioning correctly by comparing the count rate for the HAT-P-11 source spectrum and to the countrate for the binned spectrum (i.e. the “synspec” spectrum).
    The count rate for the binned spectrum is only %0.06 greater than the source spectrum’s, so I think it’s safe to say that synphot‘s binned photometry is behaving how we want it to.
  • Wrote a rough draft of what could be the synspec aspect of a higher function (see “Longer term goals” below)
I also took my dog to the groomers and I haven’t stopped squealing since

This week’s goals:

  • (minor) Fix units in SNR functions
  • Spend time with the notebook tutorials I’ve written and think about any higher order functions that could make their functionalities more streamlined
  • Take a look at telescopy and astroplan and think how they, along with synspec and any higher order functions we come up with, can be merged

Longer term goals:

  • Some thoughts on a possible deliverable, copied from my Slack message to the team:

So I’ve been thinking a bit more about how to incorporate synthetic spectroscopy into a package for synthetic photometry, and I’ve come full circle back to “telescopy”. At our last meeting we talked about writing a kind of wrapper for the synphot concepts covered in the tutorials, and I’m wondering if telescopy can be that wrapper.

Where synphot has a wide range of functionality, perhaps telescopy could be a package which utilizes synphot for a more “streamlined” purpose, for example “telescopy: a package for estimating the exposure time needed to obtain x-SNR with y-instrument”, where y-instrument could be either photometric or spectroscopic in nature.

The problem of this which we’ve discussed before is that we would be creating “yet another package”. But it would be a solution which provides both a wrapper and a synspec resource while avoiding the whole “why is there a synthetic spectroscopy generator inside of a synthetic photometry package??” confusion.

Week 8: More synspec and SNR

What I completed this week:

  • Made a pull request to synphot for implementing the signal-to-noise calculators that I wrote last week. I forgot to mention in my blog post last week that writing this module was a really useful exercise for me because I was introduced to so many new things! Like:
    • writing tests for pytest (and using the assert_quantity_allclose function from astropy.tests.helper)
    • python decorators (and how to ensure parameters are given in compatible units with astropy.units.quantity_input)
    • “duck typing” and the hasatrr() function
    • private functions
    • citing references in the docstring in the “astropythonic” way
    • a little bit of Sphinx documentation syntax
Me every day as a GSoC student
  • Wrote the first-draft of a “synspec” (synthetic spectroscopy) example. To simulate spectroscopy, I essentially create a custom waveset derived from some spectral resolution R and use synphot to bin the flux accordingly. There are a couple of things I’m unsure about though:
    • Does the binning done by synphot do (for all intents and purposes) the same thing as a real spectrometer? I haven’t done much work with spectroscopy so the mechanics are a little mysterious to me…
    • The spectral resolution R is defined to be R = λ / Δ λ , but what is the “correct” way of constructing an array of wavelengths according to this resolution? Do we keep R constant over the range of wavelengths (in which case Δ λ changes for each wavelength), or do we keep Δ λ constant and have a range of R over the wavelength range? For example, if we have an instrument which ranges from 0.3 – 1.1 microns with a spectral resolution of R = 1000, then either:
      • We let Δ λ change from 0.0003 to 0.0011 microns over the wavelength range
      • We define R = 1000 to be the resolution at the center of the waveset, (0.7 microns), which gives Δ λ = 0.0007 microns. Then for a constant Δ λ, the range of R is [430, 1600] (how it’s done in the current version of synspec)
  • Started working on the astroquery pull request for querying the SVO filter service, borrowing heavily from tynt (no need to reinvent the wheel).

This week’s goals:

  • Finish the synphot SNR pull request by adding to the docs
  • Adding SNR functionality into the other examples
  • Continue working on synspec
    • To check that the binned photometry (aka spectroscopy) is functioning correctly, we can compute the counts for a source spectrum and compare to the counts for the spectroscopized(?) spectrum. Then if they are the same within some error we can rest easier.
    • As it is now, we set a constant Δ λ and allow R to vary accordingly. Add functionality such that the user can either set a constant Δ λ and not worry about what R is, or they can give a specific R and not worry about what Δ λ is.
  • Submit the astroquery PR once it has been reviewed by mentors

Longer term goals:

  • Make 4-5 notebooks which explore different use cases in order to get an idea of how we want to implement any changes or enhancements to synphot:
    1. Model spectrum -> Synthetic photometry
      • Ground-based example: existing APO notebook
      • Space-based example: existing APO notebook + Kepler
    2. Empirical spectrum -> Synthetic photometry
      • Example: Galaxy with high H-alpha emission observed by SDSS
    3. Model spectrum -> synthetic spectroscopy
      • Example: Observations of a G dwarf with CHEOPS at R~1k vs 100k – what count rates do you get?
    4. Empirical spectrum -> synthetic spectroscopy
      • Possible example: Some HII region spectrum -> how many hours to a S/N of X
  • Make a pull request to synphot for the signal-to-noise predictor
  • Create a pull request to Astroquery to enable queries to the filter VO service for filter transmittance curves
  • Think about whether we want to add a higher-level function which takes the telescope params and target and spits out a SNR

Weeks 6 and 7: Synspec and SNR

What I completed these last two weeks:

  • Trying to solve the offset problem for example 2:
    • Do the same analysis for multiple objects. Is the systematic error still there?
      The systematic error persists strongly for both objects, as you can see in the Figure 1 below. So it’s not object-specific, and I think this also means it’s not weather related (assuming the two objects were observed on different nights).
    • Plot bandpasses on top of spectrum
      Figure 2 shows spectra of both objects ontop of the bandpasses, where we’ve removed the u ad z bands because they don’t fully overlap with the spectrum, which resulted in underestimated magnitudes for the u and z bands since synphot has to guess what the spectrum looks like to fill in the gaps. Other than that, I can’t see a reason why the filter models would be responsible for the systematic offset.
    • How systematic is this offset? Plot the fractional error in each band to see.
      Figure 3 shows the fractional error between the synphot fluxes and the observed fluxes. The offset isn’t quite as systematic as I thought, since the three objects are affected to different degrees. Could this mean that the error *is* weather related…?
    • Add a “boring” elliptical galaxy for comparison.
      The offset persists for the elliptical galaxy, and is in the same ballpark as the others in terms of fractional error. I’m still unsure why synphot would be overestimating the photometric flux…
Figure 1. Photometry of two different H-alpha emission galaxies and one elliptical galaxy observed by SDSS, with a persistant systematic offset between synphot fluxes and observed fluxes.
Figure 2. Spectra of three objects in relation to the bandpass models.
Figure 3. Fractional error between the synphot fluxes and the observed fluxes in each band.
  • Wrote a signal-to-noise predictor based on “the CCD equation” in Steve Howell’s “Handbook of CCD Astronomy”, 2000. I would like to make a pull request to implement this into synphot, possibly as a method on Observation, e.g.:
>>> from synphot.observation import Observation
>>> import astropy.units as u
...
>>> observation = Observation(spectrum, bandpass)
>>> countrate = observation.countrate()
>>> counts = countrate * 10 * u.s
>>> snr = observation.snr(counts)
>>> snr
42
  • Also wrote a function which returns the exposure time needed to reach a desired signal to noise ratio.
    • I’m not really sure what to name it, exposure_time_from_snr() seems too long…
    • In the case where the error due to the background estimation and gain are significant, solving for t can not be done analytically, so I use scipy.optimize.fsolve to solve it numerically. However, fsolve doesn’t preserve quantity types, so I have to add some annoying conditional statements to convert t back to a quantity… is there a better way to do this?

This week’s goals:

  • Make a pull request to synphot for the signal-to-noise predictor
  • Begin working on “synspec” (synthetic spectroscopy)
    • Is there an easy way to use synphot for synthetic spectroscopy?
    • See how long each run takes – ideally it should be fast, for people who want to run this on thousands of spectra

Longer term goals:

  • Make 4-5 notebooks which explore different use cases in order to get an idea of how we want to implement any changes or enhancements to synphot:
    1. Model spectrum -> Synthetic photometry
      • Ground-based example: existing APO notebook
      • Space-based example: existing APO notebook + Kepler
    2. Empirical spectrum -> Synthetic photometry
      • Example: Galaxy with high H-alpha emission observed by SDSS
    3. Model spectrum -> synthetic spectroscopy
      • Example: Observations of a G dwarf with CHEOPS at R~1k vs 100k – what count rates do you get?
    4. Empirical spectrum -> synthetic spectroscopy
      • Possible example: Some HII region spectrum -> how many hours to a S/N of X
  • Make a pull request to synphot for the signal-to-noise predictor
  • Create a pull request to Astroquery to enable queries to the filter VO service for filter transmittance curves Is this still on the table given tynt‘s existence?
Oh and it was also the 4th of July sometime in between the coding

Week 5: This time with galaxies

What I completed this last week:

  • Nearly finished the empirical spectrum -> synthetic photometry tutorial for a high H-alpha emission galaxy. This example can be seen here.
    • As you can see from the 1-1 plot at the end of the tutorial, the predicted photometry does not perform as well compared to the 1-1 plot for the first use-case. As Erik (one of my mentors, for anyone reading this who is not already one of my mentors) mentioned, even having a percentage error within 30% is not terrible if you just want a ballpark estimate of the required exposure time.
    • This systematic offset still bugs me though, because it seems like it could be easily solved:
I suppose the mystery continues, for now…
  • Took part in an astropy hack day at the Center for Computational Astrophysics in NYC! Erik and I got to work together in person, I met some astronomers I knew of but never actually met, and I got to witness a bit more of the inner workings of open source development. Also, free food. Lots of free food.

This week’s goals:

  • Try to solve the offset problem for example 2
    • Do the same analysis for multiple objects. Is the systematic error still there?
    • Plot bandpasses on top of spectrum
  • Begin working on synthetic *spectroscopy* (synspec)?
    • Is there an easy way to use synphot for synthetic spectroscopy?
    • See how long each run takes – ideally it should be fast, for people who want to run this on thousands of spectra
  • Begin working on a signal-to-noise predictor (maybe as an eventual method on Observation?)
Me when there’s free food

Longer term goals:

  • Make 4-5 notebooks which explore different use cases in order to get an idea of how we want to implement any changes or enhancements to synphot:
    1. Model spectrum -> Synthetic photometry
      • Ground-based example: existing APO notebook
      • Space-based example: existing APO notebook + Kepler
    2. Empirical spectrum -> Synthetic photometry
      • Example: Galaxy with high H-alpha emission observed by SDSS
    3. Model spectrum -> synthetic spectroscopy
      • Example: Observations of a G dwarf with CHEOPS at R~1k vs 100k – what count rates do you get?
    4. Empirical spectrum -> synthetic spectroscopy
      • Possible example: Some HII region spectrum -> how many hours to a S/N of X
  • Implement a signal-to-noise predictor
  • Create a pull request to Astroquery to enable queries to the filter VO service for filter transmittance curves

Week 4: SkyCalc and SVO and tynt oh my

What I completed this last week:

  • Successfully queried from the SVO filter database for our tutorial’s throughput models. This inspired Brett to write a package called tynt, a “super lightweight package containing approximate transmittance curves for more than five hundred astronomical filters”, which I’ve implemented into the synphot examples.
  • Added a Kepler (i.e. space-based) example to our tutorial. The errors between the synphot counts and the empirical counts are less than 15%.
  • Wrote atmospheric_transmittance.py to model the atmospheric transmittance with the SkyCalc Sky Model Calculator. Parameters can be set when calling get() for a more precise model of the sky, otherwise get() will use the default parameters provided by SkyCalc.

This week’s goals:

  • Begin working on Example 2: Empirical spectrum (like from SDSS/Hubble website) -> Synthetic photometry
  • Begin working on the Astroquery pull request? (mentioned below)

Longer term goals:

  • Make 4-5 notebooks which explore different use cases in order to get an idea of how we want to implement any changes or enhancements to synphot:
    1. Model spectrum -> Synthetic photometry
      • Ground-based example: existing APO notebook
      • Space-based example: existing APO notebook + Kepler
    2. Empirical spectrum (like from SDSS/Hubble website) -> Synthetic photometry
      • Example: Erik’s palomar spectrum + MDM Halpha observations
    3. Model spectrum -> synthetic spectroscopy
      • Example: Observations of a G dwarf with [the space-based mission Brett mentioned called CHEOPS, but don’t worry about that] at R~1k vs 100k – what count rates do you get?
    4. Empirical spectrum -> synthetic spectroscopy
      • Possible example: Some HII region spectrum -> how many hours to a S/N of X
  • Implement a signal-to-noise predictor
  • Create a pull request to Astroquery to enable queries to the filter VO service for filter transmittance curves

Week 3: Querying progress

What I completed last week:

  • Made a separate git repo for my GSoC related work. This will make it easier to track changes, make suggestions, etc
  • Edited the first tutorial such that:
    • The only data files left in notebook are the CCD quantum efficiency table. The other bits of data are queried/modeled with astropy.utils.download_file(), astroquery.gaia, and pwv_kpno.
    • There is now a preamble including authors, objectives, keywords, and links to different sections throughout the tutorial
  • Contacted the Spanish Virtual Observatory and asked them to add the primed SDSS filters (which are used at APO) to their filter profile service, which they have now done!
  • Explored using pwv_kpno as an atmospheric model. I was struggling with it at first, but I didn’t wait too long to raise an issue on its github page. Daniel was very helpful and responsive, and things are running smoothly now. Because pwv_kpno‘s main functionality is modeling the effects of precipitable water vapor on the atmospheric transmission, it isn’t a complete atmospheric model (it doesn’t claim to be) since it doesn’t include the opacity due to Rayleigh scattering. While this would be okay for observations in the infrared, it significantly affects our visual-band count estimates.
  • Addressing the above, I’m going to try Brett’s suggestion to use skycalc_cli and see if that makes a more complete atmospheric model. To do this, I am going to borrow the bit about querying from Cerro Paranal from skycalc, and perhaps eventully turn this into an astroquery pull request.

This week’s goals:

  • In the first tutorial, edit the bandpass retrieval to query from SVO instead of APO. Alert the github universe when it’s ready to be looked over!
  • Create a short example of Kepler (i.e. space-based) counts for HAT-P-11 and TRAPPIST-1
  • Keep investigating skycalc_cli for atmospheric transmission models
  • Begin working on Example 2: Empirical spectrum (like from SDSS/Hubble website) -> Synthetic photometryExample: Erik’s palomar spectrum + MDM Halpha observations
  • Begin working on the Astroquery pull request? (mentioned below)

Longer term goals:

  • Make 4-5 notebooks which explore different use cases in order to get an idea of how we want to implement any changes or enhancements to synphot:
    1. Model spectrum -> Synthetic photometry
      • Ground-based example: existing APO notebook
      • Space-based example: existing APO notebook + Kepler
    2. Empirical spectrum (like from SDSS/Hubble website) -> Synthetic photometry
      • Example: Erik’s palomar spectrum + MDM Halpha observations
    3. Model spectrum -> synthetic spectroscopy
      • Example: Observations of a G dwarf with [the space-based mission Brett mentioned called CHEOPS, but don’t worry about that] at R~1k vs 100k – what count rates do you get?
    4. Empirical spectrum -> synthetic spectroscopy
      • Possible example: Some HII region spectrum -> how many hours to a S/N of X
  • Implement a signal-to-noise predictor
  • Create a pull request to Astroquery to enable queries to the filter VO service for filter transmittance curves
That’s a lot of water vapor…

Week 2: Updates and a second to-do list

We’ve made some headway with the synphot tutorial, which has helped us determine what to do next. This is what we’ve done so far:

  • To get more accurate count rates, we:
    • Add effects by atmospheric attentuation by using the Cerro Paranal model transmittance curves for an airmass of 1.5
    • Consider the effects of the quantum efficiency on the spectra by using the values in the table found in section 3.5 on this page of APO’s website
    • Model the source spectra using model spectra from PHOENIX instead of blackbody models. For HAT-P-11 we use Teff = 4800 K, and for TRAPPIST-1 we use Teff = 2500 K , both with logg = 4.5 cm / s^2
    • Realized that we have to divide the output of synphot’s countrate() function by the gain of the modeled telescope
  • There was some debugging we had to tackle to obtain the correct units/order of magnitude for the PHOENIX source spectra – the blackbody model is in units of 𝑒𝑟𝑔 𝑠−1 𝑐𝑚−2 𝐴˚−1 𝑠𝑟−1, while PHOENIX gives flux in 𝑒𝑟𝑔 𝑠−1 𝑐𝑚−3. I don’t think I fully understand yet, but the problem seems to be that while synphot handled the cm to angstrom conversion fine, the steradian was sort of lost in translation… The SourceSpectrum object was correct as long as we divided the normalized flux by pi.

With these factors implemented, our current precision looks like:

synphot TRAPPIST-1: 257K
actual TRAPPIST-1: 203K
synphot HAT-P-11: 30M
actual HAT-P-11: 34M

This week’s goals:

  • Make a separate repo for notebook tutorials (maybe make a github “project” out of it?)
  • Have no data files left in notebook (except CCD QE), get the needed data by querying instead
    • In the meantime we will use this new functionality to query the SDSS filters for our count rate example
  • Investigate using pwv_kpno to compute transmittance rather than using the Cerro Paranal model to further improve count rates.

Longer term goals:

  • Make 4-5 notebooks which explore different use cases in order to get an idea of how we want to implement any changes or enhancements to synphot:
    1. Model spectrum -> Synthetic photometry
      • Example: existing APO notebook
      • (And maybe also Kepler?)
    2. Empirical spectrum (like from SDSS/Hubble website) -> Synthetic photometry
      • Example: Erik’s palomar spectrum + MDM Halpha observations
    3. Model spectrum -> synthetic spectroscopy
      • Example: Observations of a G dwarf with [the space-based mission Brett mentioned called CHEOPS, but don’t worry about that] at R~1k vs 100k – what count rates do you get?
    4. Empirical spectrum -> synthetic spectroscopy
      • Possible example: Some HII region spectrum -> how many hours to a S/N of X
  • Implement a signal-to-noise predictor