Flux bug fix in quicktransients simulator

[sybenzvi]

Fixed bug where a small number of spectra generated with quicktransients had zero flux. It was due to a mistake in quicktransients.py; an incorrect redshift argument was being passed to sim_spectra. Fixes issue #540.

While poking around in quicktransients.py, did some additional cleanup:

• argparse help now prints defaults
• use the logger instead of print statements
• coadd camera spectra if requested

[sybenzvi]

Note: the failing check in travis seems related to astropy.

[sbailey]

@moustakas could you vet the usage of redshift here? If I understand correctly, this change means that the truth table redshifts are never used for quicktransients, which isn't obviously right to me (vs. fixing why using those redshifts for BGS hosts was resulting in zero flux spectra).

I'm checking on the Travis astropy failures; we're pinning the astropy version to 2.0.16 (or at least trying to...) and that version with this branch passes at NERSC, so it isn't yet clear to me what changed at Travis.

[moustakas]

I'll review when I can but I'm focused on DR9 right now.

[sybenzvi]

@sbailey, I'm certain the truth table redshifts are used correctly in the transients simulation because the generated template spectra are sampled and added in the host rest frame and then the summed template is returned in the observer frame.

The redshift argument described here is an array passed to desisim.simexp.simulate_spectra. The values are used to adjust the half-light radii of targets during the simulation (line 558 of desisim.simexp.py) of BGS spectra. That's why they're not applied to other galaxy types.

It seems like a few percent of the time the size adjustment produces a divide-by-zero during the spectrum generation. I didn't understand why and instead opted to set redshift=None in quicktransients.py... but that was lazy on my part. Since @moustakas is preoccupied I'll dig a little and see if I can track down where it goes wrong.

[sybenzvi]

I looked into it and what's happening is:

1. If redshift is not Nonein the argument ofdesisim.simexp.simulate_spectra` it is used to adjust the angular size of BGS galaxy bulges and disks (line 560 of desisim/simexp.py). When the simulated redshift is small, e.g., z < 0.01, the size gets quite large.
2. In specsim, the resulting fiber acceptance fraction becomes zero for these cases.
3. The source_flux and source_flux_to_photons is then set to zero (line 438, 454 of specim/simulator.py).
4. The output flux_calibration (line 549) divides by source_flux_to_photons and emits a divide by zero runtime warning.

The quicktransients script is using the default instrument settings in specsim, which means it's using fastsim. Maybe that is not compatible with the redshift argument?

I noticed that quickspectra does not pass a redshift array to desisim.simexp.simulate_spectra, so maybe this PR is correct to also remove it from quicktransients. What do you think @sbailey?

[sbailey]

Thank for digging into this @sybenzvi . It looks like the fundamental issue is that specsim fastsim allows the fiber acceptance fraction to go to zero instead of just a very small number. This could be fixed in specsim, or worked around in desisim by capping the size of the galaxies it considers. I'm not sure if it properly tracks the total flux for galaxies that become gigantinormous either.

It should be more correct for the host galaxy sims to use the truth redshift than for them to not use it, and I think desisim.simexp.simulate_spectra not using the redshift array was a shortcut rather than correct behavior to be followed. The redshift and fiberlosses matter more for big BGS galaxies than small ELG/LRG, which is a partial explanation for simulate_spectra not using the redshifts (still not right, but less critical).

FWIW, the travis problem was that recent healpy is incompatible with old astropy=2.x, but its dependencies don't list that correctly so Travis was installing an version of healpy that was incompatible to our pinned astropy=2.x. We're working on upgrading our full stack to astropy=4.x, but don't have the full suite working yet to then upgrade the tests too (goal is to go through one set of tags+software release that supports both astropy=2 py=3.6 and astropy=4 py=3.8 and then deprecate support for the old versions).

[sybenzvi]

@sbailey, I realized that I never checked if the zero flux issue also depends on observing conditions. I always simulated objects with 1.1 arcsec seeing. Before committing to one solution we should cap the angular size of objects in desisim.simexp.simulate_spectra and vary the seeing. If the zero flux issue becomes more common as the seeing gets worse, then we know the problem is in specsim and should be addressed there. If the angular size cap works regardless of simulated observing conditions, then we can cap object sizes and push the change to this branch of desisim. I'll be preoccupied through tomorrow afternoon but should have it done by Thursday.

[moustakas]

This PR has been hanging long enough that I'd like to get it merged. Let's continue discussion of any outstanding issues in a new ticket. (I'm also concerned that #559 may have broken some of the transient simulation code, which I didn't get a chance to investigate in that PR; sorry @sybenzvi!)

[coveralls]

Coverage decreased (-0.08%) to 45.628% when pulling 20c07c0 on fluxbug_fix into 9440c3c on master.