Implement more flexible fiber acceptance fraction calculations

[dkirkby]

We currently use tabulated fiber acceptance fractions for a fixed set of assumed transverse source profiles on the sky. This feature request is to allow arbitrary source profiles to be simulated, probably parameterized by a sum of Sersic components.

[moustakas]

Stephen and I have had some discussions about this. The spectral templates I've developed don't yet--but will---have some morphological information attached to them (half-light radii, Sersic indices, axis ratios, etc.), at least for ELGs, LRGs, and BGS (based on HST observations wherever possible). The N-dimensional space of size, SB profile, axis ratio, fiber size (which depends on focal-plane position), seeing, fiber mis-centering, etc. is quite large, however, so this will require some care as to how these effects will be implemented. +1 that this is being tackled and I'm happy to help!

[dkirkby]

Most of the work is already done in this exposure-time calculator notebook, so this is mostly a matter of repackaging existing code. The profile and PSF modeling are handled by GalSim, so are very flexible. However, I don't want to make GalSim a required dependency so we will continue to support the use of tabulated fiber acceptance ratios for predetermined profiles, and only use GalSim for more flexible calculations when it is available.

[sbailey]

Quick notes:

This issue and also #3 (guiding and tracking errors) are things we'd like to simulate for the input to the full pixsims as well. Please implement these in a way that they can be used independent from the rest of the specsim simulation, i.e. as separate library functions and not embedded within a larger calculation.

DESI has a radially dependent platescale, which results in elliptical fiber profiles in angular coordinates on the sky. Ideally this code will support not only arbitrary source profiles, but also non-circular apertures.

For reference, desimodel/bin/fiberloss.py is the code we used to calculate the current fiber input geometric loss factors. It brute forces the calculation over a fine grid, which is viable if you are doing it once but not fast enough if you are doing it for every fiber on the fly for a custom profile.

[dkirkby]

We are already modeling the SDSS radial distortions for the spectrophoto calibration paper, so that sounds doable. Is the DESI radial distortion function already tabulated somewhere in desimodel?

[dkirkby]

Please implement these in a way that they can be used independent from the rest of the specsim
simulation, i.e. as separate library functions and not embedded within a larger calculation.

That's definitely the plan. Another piece that should be easy to call independently is converting flux and throughput arrays to an array of mean photon rates, which I see is already being duplicated in quicksim and specter.throughput.Througput.photons() (which desisim.obs.new_exposure calls via desimodel.io.load_throughput, creating some interesting dependencies!)

[sbailey]

The DESI optics platescale vs. radius is tabulated in data/focalplane/platescale.txt, coming from DESI-0329v14 which has more explanation. Note that the latest version is v15, but the platescale info is the same.

Throughput calculation: in the short term it may be useful for us to have independent paths to compare and debug. After specsim matures, I could make this a dependency of specter depending upon how many things like throughput are an overlap of functionality.

[moustakas]

Any plans for moving the resolution of this issue up? I'm working on simulating BGS objects (as part of desihub/desisim#72) but specsim only knows about object types 'lrg', 'elg', 'qso' and 'star'.

Or should I just make a nominal 'fiberloss-bgs.dat' file? The problem with having just one file, of course, is that we expect significantly larger variations in the fiber acceptance rate for BGS because of the flux-limited and lower-redshift nature of the sample. It would be great for my simulations to marginalize over these effects.

[dkirkby]

I tagged this issue with the v0.5 milestone, which means that I plan to address it in the next release.

Once this is implemented, you will need some model of BGS shape profiles (e.g., in half-light radius and Sersic index) in order to sample the variability in fiber loss fractions.

[dkirkby]

@moustakas You wrote earlier:

...we expect significantly larger variations in the fiber acceptance rate for BGS because
of the flux-limited and lower-redshift nature of the sample.

I am going to model the dependence of fiber acceptance fraction on the source profile. Will that capture the BGS variations you expect? How flexible a profile model do you need?

• Can I assume that all profiles are circularly symmetric, or will source ellipticity be an important consideration?
• Is a single-component Sersic model sufficient to capture variations in BGS profiles? Can we fix Sersic n=1 so the only shape parameter is the half-light radius?

[dkirkby]

The original fiberloss calculator makes the following assumptions:

# Fiber diameter 1.5 arcsec
# Seeing is 1.1 arcsec FWHM with Moffat beta=3.5, and scale with
# wavelength as (6355./wave)^0.2
# Lateral offset is 0.2 arcsec from fiber center
# Cases are:
#   * Point source
#   * ELG with half-light radius 0.45 arcsec
#   * LRG with half-light radius 1.0 arcsec

The resulting curves are:

[moustakas]

* Can I assume that all profiles are circularly symmetric, or will source ellipticity be 
an important consideration?
* Is a single-component Sersic model sufficient to capture variations in BGS profiles? 
Can we fix Sersic n=1 so the only shape parameter is the half-light radius?

Ellipticity, the Sersic indices of the bulge and disk, half-light radii, bulge-to-disk ratio, etc., will all be important for the fiber acceptance fraction (and not just the scatter -- the mean fraction will vary in important ways). There's no way -- for BGS at least -- a single n=1 profile of varying half-light radius will be flexible enough. (Remember that the BGS is flux-limited, so we'll include both star-forming disk and spheroidal galaxies in the sample.) The question is how many effects do we need to model to capture the mean trends (and variance) of the galaxy population. There's ton's of data we can use (from SDSS, AGES, GAMA, etc.), but let's make a plan on-list.

[dkirkby]

The canonical galaxy sources are modeled as:

• ELG: n = 1 Sersic with hlr = 0.45"
• LRG: n = 4 Sersic with hlr = 1.0"
  Specifically:

img_elg    = numpy.exp(-1.678*radius/r_elg)
img_lrg = numpy.exp(-7.67*((radius/r_lrg)**0.25 - 1))

Note that both profiles are assumed to be circular, which is more of an edge case than typical.

[dkirkby]

Seeing is fixed at FWHM=1.1" at 6355A, which is consistent with Arjun & Dey 2013, specifically these panels of their Fig.11:

These seeing values include an instrumental PSF specified in desi.yaml as:

#- Jacoby seeing (subtract in quadrature from Dey and Valdes seeing values)
jacoby_seeing: 0.219         # arcsec FWHM

and subtracted in quadrature in fiberloss.py:

seeing_sigma = numpy.sqrt( (seeing/2.35482)**2 - mayall_blur**2 )
site_seeing = seeing_sigma * 2.35482

The atmospheric PSF is modeled as a Moffat profile with beta=3.5 and a FWHM that scales with wavelength with the Kolmogorov -1/5 exponent. Perhaps better to use a Kolmogorov PSF directly?

[dkirkby]

The instrumental contribution to the PSF consists of an optical PSF and distribution of expected centroid offsets, tabulated in $DESIMODEL/data/inputs/throughput/DESI-0347-throughput.txt:

# Lateral offset is sigma microns                               
# Blur is 1-d sigma microns                             

This is implemented in fiberloss.py as a Gaussian PSF with:

sigma=blur/platescale
platescale= 107.0 um / 1.52" ~ 70.4 um/arcsec

Note that this neglects the variation of platescale over the field of view. Is a Gaussian PSF the appropriate model or should it be more like the Airy pattern of an obscured pupil? What is responsible for the wavelength dependence here? (The Airy PSF size increases linearly with wavelength but is about an order of magnitude smaller than these values).

The offset is fixed at 0.2" ~ 14.1 um, which seems high given the tabulated RMS values ~ 10um. Perhaps it would be better to sample a distribution of offsets via an additional PSF convolution? It would also be feasible to simulate offsets during the exposure, a la Margala 2015, assuming they are dominated by residual ADR (with a fixed ADC) rather than guiding errors.

[dkirkby]

The DESI-0347 spreadsheet contains a plot with optical PSF size vs field angle at different wavelengths (from zemax?):

The variation with field angle is presumably in addition to the plate scale variation. I guess the desimodel blur values are some sort of average over the field of view (area weighted?).

The lower plot shows the optical dispersion over the field of view and is a potentially large effect on the fiber acceptance fraction. How is this translated to an RMS offset value in desimodel?

[djschlegel]

The fibers + spectrographs should accept almost all of the light from infinity, except for fibers that have non-compliant FRD (of which I'm sure there'll be some). The PSFs on the focal plane were generated using a N thousand ray traces through Zemax, and then the RMS is computed from the results.

The plate scale variations are large. There's both a 10% change in the projected area from the center of the focal plane to the edge, and a squashing of the projected area on the sky by 10%. Those data here: https://desi.lbl.gov/trac/browser/code/desimodel/trunk/data/focalplane/platescale.txt
ces of
It would be excellent to think about how one would recover the PSF as seen on-sky for each object, which could then be convolved with the galaxy images in the case of extended sources. The relevant pieces of information are:
(1) PSF from the guider images (backing out the squashing from these being at the edge of field), retaining any guider offsets from the commanded star locations (otherwise track separately)
(2) Chromatic distortions, which should mostly be a broadening of the PSF as one goes blue; the guider is only r-band
(3) Optical distortion at the location of each fiber
(4) convolve with the Tractor model parameters of galaxies

[dkirkby]

Is the Tractor model of the galaxy profiles documented somewhere?

[djschlegel]

http://legacysurvey.org/dr2/catalogs/

[dkirkby]

Has someone already made plots of the Tractor EXP+DEV shape parameters for the different DESI target classes?

[dkirkby]

It looks like the original fiberloss.py is missing a fwhm/sigma ~ 2.35 conversion factor here:

seeing_sigma = numpy.sqrt( (seeing/2.35482)**2 - mayall_blur**2 )

or else this comment in desi.yaml identifying the Jacoby seeing as a FHMW value is incorrect:

jacoby_seeing: 0.219         # arcsec FWHM

With the default 1.1" seeing, this is a ~10% effect on the size of the atmospheric PSF.

[sbailey]

For clarification, the mayall_blur variable in the code is read from the jacoby_seeing value in desi.yaml. Digging backwards, this is more confusing than I thought it would be: That value is derived in more detail on a tab of the spreadsheet in DESI-0347 which also refers to it as FWHM. This is based on figure 4 of Jacoby et al where it is called "RMS full-width", which seems to exactly split the confusion between whether this is an RMS (i.e. sigma-like) or a Full Width Half Max. Or maybe "RMS full-width" means 2sigma?

Comparing to the DESI corrector blur numbers calculated as a sigma (not FWHM), the numbers are about the same (not a factor of ~2 different), which under the assumption that the DESI corrector blur is not radically better/worse than the current Jacoby corrector, would imply that the Jacoby blur numbers are actually sigmas. I chatted with Pat Jelinsky and he may have the Zemax model of the current corrector, in which case he can try to reproduce figure 4 to figure out the units...

[moustakas]

Has someone already made plots of the Tractor EXP+DEV shape parameters for the different 
DESI target classes?

To my knowledge this has not been done. Briefly, the best ELG truth region is in the north (DEEP2/Field 1), and therefore does not have DECaLS data, and the LRG target selection has been somewhat of a moving target (pun intended). The morphologies and shapes of a BGS sample could be, but have not been looked at. That said we do have HST morphologies (Sersic indices, axis ratios, and half-light radii) of a few hundred ELG- and LRG-like targets that could be useful, and for BGS we basically have a ton of data.

Should I task someone on the Targeting WG to put these pieces together?

[dkirkby]

Closing since the original issue (more flexible fiberloss calculations) has been resolved since version 0.8.