desisim.simexp.simulate_spectra not varying exposure time

[humnaawan]

Hi, I have been trying to figure out why my spectra don't change with exposure time, and I've isolated the issue to desisim.simexp.simulate_spectra.

Specifically, say I have a wave array wave_arr and flux array flux_arr, and I run

simulator = desisim.simexp.simulate_spectra(wave=wave_arr,
                                            flux=flux_arr,
                                           )
out0 = simulator.camera_output[0]

vs.

simulator = desisim.simexp.simulate_spectra(wave=wave_arr,
                                            flux=flux_arr,
                                            obsconditions=obsconditions_dict
                                           )

out1 = simulator.camera_output[0]

where

obsconditions_dict = 
{'TILEID': 11147,
 'RA': 156.5882658712762,
 'DEC': 22.643639342687276,
 'PASS': 1,
 'EXPID': 344,
 'MJD': 58824.50497106333,
 'EXPTIME': 4700,
 'SNR2FRAC': 0.5060822,
 'AIRMASS': 1.0281637,
 'SEEING': 0.7948105,
 'TRANSP': 0.82818073,
 'SKY': 1.0,
 'PROGRAM': 'DARK',
 'NIGHT': '20191206',
 'FLAVOR': 'science',
 'MOONFRAC': 0.5,
 'MOONALT': -10.0,
 'MOONSEP': 90.0}

Then I'd expect out0 and out1 to not be the same -- is that not right?

The calls to simexp seems to suggest that the observation parameters being considered are indeed the ones I pass, as I see the following for the first run

DEBUG:simexp.py:418:simulate_spectra: loading specsim desi config desi
DEBUG:simexp.py:422:simulate_spectra: creating specsim desi simulator
WARNING:simexp.py:428:simulate_spectra: Assuming DARK conditions
INFO:simexp.py:444:simulate_spectra: MJD not in obsconditions, using DATE-OBS 2019-12-07T12:07:09.500
DEBUG:simexp.py:448:simulate_spectra: obsconditions SEEING = 1.1
DEBUG:simexp.py:448:simulate_spectra: obsconditions EXPTIME = 1000
DEBUG:simexp.py:448:simulate_spectra: obsconditions AIRMASS = 1.0
DEBUG:simexp.py:448:simulate_spectra: obsconditions MOONFRAC = 0.0
DEBUG:simexp.py:448:simulate_spectra: obsconditions MOONALT = -60
DEBUG:simexp.py:448:simulate_spectra: obsconditions MOONSEP = 180
WARNING:simexp.py:653:get_source_types: (SV1_)DESI_TARGET not in fibermap table; using source_type='star' for everything
DEBUG:simexp.py:490:simulate_spectra: running simulation with fastsim fiber loss method
DEBUG:simexp.py:496:simulate_spectra: source types: 500 star

and this for the second

DEBUG:simexp.py:418:simulate_spectra: loading specsim desi config desi
DEBUG:simexp.py:422:simulate_spectra: creating specsim desi simulator
INFO:simexp.py:442:simulate_spectra: exposure_start 2019-12-07T12:07:09.500
DEBUG:simexp.py:448:simulate_spectra: obsconditions SEEING = 0.7948104739189148
DEBUG:simexp.py:448:simulate_spectra: obsconditions EXPTIME = 4700
DEBUG:simexp.py:448:simulate_spectra: obsconditions AIRMASS = 1.0281636714935303
DEBUG:simexp.py:448:simulate_spectra: obsconditions MOONFRAC = 0.5
DEBUG:simexp.py:448:simulate_spectra: obsconditions MOONALT = -10.0
DEBUG:simexp.py:448:simulate_spectra: obsconditions MOONSEP = 90.0
WARNING:simexp.py:653:get_source_types: (SV1_)DESI_TARGET not in fibermap table; using source_type='star' for everything
DEBUG:simexp.py:490:simulate_spectra: running simulation with fastsim fiber loss method
DEBUG:simexp.py:496:simulate_spectra: source types: 500 star

but then I don't understand why the output spectra are not different. I'd really appreciate some help since I don't know what is going wrong. I thought it was my fibermap (that I was passing to desisim.simexp.simulate_spectra) but that seems to be making no difference to the actual outputs as well, although I do see that the source types are read:

DEBUG:simexp.py:418:simulate_spectra: loading specsim desi config desi
DEBUG:simexp.py:422:simulate_spectra: creating specsim desi simulator
INFO:simexp.py:442:simulate_spectra: exposure_start 2019-12-07T12:07:09.500
DEBUG:simexp.py:448:simulate_spectra: obsconditions SEEING = 0.7948104739189148
DEBUG:simexp.py:448:simulate_spectra: obsconditions EXPTIME = 4700
DEBUG:simexp.py:448:simulate_spectra: obsconditions AIRMASS = 1.0281636714935303
DEBUG:simexp.py:448:simulate_spectra: obsconditions MOONFRAC = 0.5
DEBUG:simexp.py:448:simulate_spectra: obsconditions MOONALT = -10.0
DEBUG:simexp.py:448:simulate_spectra: obsconditions MOONSEP = 90.0
DEBUG:simexp.py:490:simulate_spectra: running simulation with fastsim fiber loss method
DEBUG:simexp.py:496:simulate_spectra: source types: 180 elg 10 star 167 qso 103 lrg 40 sky
WARNING:simexp.py:511:simulate_spectra: the half light radii are fixed here for LRGs and ELGs (and not magnitude or redshift dependent)

The only other thing I can think of is the config file but then I don't understand why there's an obsconditions optional parameter at all if I am supposed to be passing a non-default config file (which would be a pain to create for every tile).

Thanks so much!

[moustakas]

Thanks for the ticket, @humnaawan, I can try to look into this.

[andreufont]

Thanks @moustakas for looking into this. At least in the Lyman alpha WG none had played with these observing conditions before, so I'm not too surprised that there are hiccupts. I hope these are easy to fix!

[alxogm]

Hi all, thanks @moustakas for taking a look.

For reference, in the Y1 simulations we make, we do vary the exposure time, and it does change the spectra. It uses the same routine desisim.simexp.simulate_spectra Humna is trying to use.

I see in quickquasars we start with a reference set of obsconditions, which we then modify, here is the line we start with,

desisim/py/desisim/scripts/quickquasars.py

Line 874 in d9b552f

obsconditions = reference_conditions[args.program.upper()]

we then just modify the exptime for each target as needed, and I think the same applies for the other arguments in obsconditions. So maybe is issue at the dictionary level?

[andreufont]

@moustakas - any update on this? We will discuss this project at the KP6 telecon later today.

[moustakas]

Apologies, I have not had a chance to look yet. But @humnaawan can you include an end-to-end example of what you're doing, including the input spectra you're passing to simulate_spectra? Thanks!

[humnaawan]

no worries @moustakas! here's the code that should produce the issue:

import numpy as np
from astropy.io import fits
from astropy.table import Table
import desisim.simexp

# read in the simpec file: wave, flux, obsconditions
simspec_fname = '/global/cfs/cdirs/desi/users/awan/20191206/00000344/simspec-00000344.fits'
simspec_hdul = fits.open(simspec_fname)
wave_simspec = simspec_hdul['WAVE'].data
flux_all_fibers = simspec_hdul['FLUX'].data
obsconditions = simspec_hdul['OBSCONDITIONS'].data
obsconditions_dict = [dict(zip(obsconditions.columns.names, row)) for row in obsconditions][0]

# read in the fibermap
fibermap_fname = '/global/cfs/cdirs/desi/users/awan/20191206/00000344/fibermap-00000344.fits'
fibermap = Table.read(fibermap_fname, format='fits')
fibermap_by_petal = fibermap.group_by('PETAL_LOC')

# get the flux for one petal
petal_ind = 0
fibermap_this_petal = fibermap_by_petal.groups[petal_ind]
flux_this_petal = flux_all_fibers[fibermap_this_petal['FIBER'].quantity, :] # assuming row index = fiber number

# set up the simulator object without any observing conditions
simulator = desisim.simexp.simulate_spectra(wave=wave_simspec,
                                            flux=flux_this_petal,
                                           )
out0 = simulator.camera_output

# now set up the simulator object with specific observing conditions
simulator = desisim.simexp.simulate_spectra(wave=wave_simspec,
                                            flux=flux_this_petal,
                                            obsconditions=obsconditions_dict
                                           )

out1 = simulator.camera_output

# compare the outputs
for cam in range(3):
    print(f'things match: {np.unique(out0[cam] == out1[cam])}')

in case its helpful, i've put a copy of the notebook with all the code above at /global/cfs/cdirs/desi/users/awan/notebooks/. please lmk if this code is not able to reproduce the issue with you (i.e. out0 is exactly the same as out1).

[sybenzvi]

@humnaawan, I get the same results as you with the current main branches of desisim and specsim.

However, I decided to see if the full sim_spectra function in desisim.scripts.quickspectra responds correctly to changes in EXPTIME or is insensitive as you reported before. For me the script behaves reasonably. I increased the EXPTIME to an absurd number to make differences in the fluxes obvious, and I do see a difference. I'm loading the desiconda master package when running this.

from desisim.scripts.quickspectra import sim_spectra
from desispec.io import read_spectra
...

obsconditions_dict['EXPTIME'] = 100000.0

# Run specsim, generate realizations of the flux model, and output to a FITS file
sim_spectra(wave=wave_simspec,
            flux=flux_this_petal, 
            program='dark', 
            obsconditions=obsconditions_dict, 
            spectra_filename='out1.fits')

spec1 = read_spectra('out1.fits')

# plot stuff...

[alxogm]

@humnaawan I think the issue is that the simulator object is not being cleared by just reusing the variable with different arguments. I added this line to your notebook

desisim.specsim._simulators.clear()
before calling the simulator for the second time and that cleared everything and returned different spectra. Can you try it also fix the issue for you?

[humnaawan]

@sybenzvi thanks! so i was largely following how spectra are created in quickspectra which is why i was stumped at what i was seeing; the reason why i'm working with the simulator object directly is that we need to handle instrument response in more detail than is done in quickspectra.

adding desisim.specsim._simulators.clear() as @alxogm suggested works though so that's great. maybe this should be mentioned somewhere since the documentation for desisim.simexp.simulate_spectra gives no impression that the simulator object needs to be cleared for things to not be wrong.