Merge pull request #2337 from desihub/variablethru

fix usage of FIBERSTATUS VARIABLETHRU

py/desispec/fiberbitmasking.py

@@ -147,7 +147,7 @@ def get_skysub_fiberbitmask_val(band):
     Return mask of bad FIBERSTATUS bits for selecting sky fibers,
     i.e. fibers with these bits set should not be used for the sky model
     """
-    return get_all_fiberbitmask_with_amp(band) | fmsk.VARIABLETHRU
+    return get_all_fiberbitmask_with_amp(band)
 
 def get_flat_fiberbitmask_val(band):
     """
@@ -169,7 +169,7 @@ def get_stdstars_fiberbitmask_val(band):
     Return mask of bad FIBERSTATUS bits for selecting standard stars,
     i.e. fibers with these bits set should not be used as standard stars
     """
-    return get_all_fiberbitmask_with_amp(band) | fmsk.POORPOSITION | fmsk.VARIABLETHRU
+    return get_all_fiberbitmask_with_amp(band) | fmsk.POORPOSITION
 
 def get_all_nonamp_fiberbitmask_val():
     """Return a mask for all fatally bad FIBERSTATUS bits except BADAMPB/R/Z

py/desispec/fluxcalibration.py

@@ -18,6 +18,7 @@
 from desitarget.targets import main_cmx_or_sv
 from desispec.fiberfluxcorr import flat_to_psf_flux_correction,psf_to_fiber_flux_correction
 from desispec.gpu import is_gpu_available, NoGPU
+from desispec.maskbits import fibermask
 import scipy, scipy.sparse, scipy.ndimage
 import sys
 import time
@@ -89,6 +90,9 @@ def isStdStar(fibermap, bright=None):
             log.warning('Using FA_TYPE to find standard stars instead')
             yes = (fibermap['FA_TYPE'] & FA_STDSTAR_MASK) != 0
 
+    #- Remove fibers with known variable throughput
+    yes &= ((fibermap['FIBERSTATUS'] & fibermask.VARIABLETHRU) == 0)
+
     return yes
 
 def applySmoothingFilter(flux,width=200) :

py/desispec/sky.py

@@ -21,7 +21,7 @@
 import sys
 from desispec.fiberbitmasking import get_fiberbitmasked_frame_arrays, get_fiberbitmasked_frame
 import scipy.ndimage
-from desispec.maskbits import specmask
+from desispec.maskbits import specmask, fibermask
 from desispec.preproc import get_amp_ids,parse_sec_keyword
 from desispec.io import findfile,read_xytraceset
 from desispec.calibfinder import CalibFinder
@@ -198,6 +198,8 @@ def get_sky_fibers(fibermap, override_sky_targetids=None, exclude_sky_targetids=
     some targetids by providing a list of them through exclude_sky_targetids
     or by just providing all the sky targetids directly (in that case
     the OBJTYPE information is ignored)
+
+    Fibers with FIBERSTATUS bit VARIABLETHRU are also excluded
     """
     log = get_logger()
     # Grab sky fibers on this frame
@@ -206,7 +208,9 @@ def get_sky_fibers(fibermap, override_sky_targetids=None, exclude_sky_targetids=
         skyfibers = np.where(np.in1d(fibermap['TARGETID'], override_sky_targetids))[0]
         # we ignore OBJTYPEs
     else:
-        skyfibers = np.where(fibermap['OBJTYPE'] == 'SKY')[0]
+        oksky = (fibermap['OBJTYPE'] == 'SKY')
+        oksky &= ((fibermap['FIBERSTATUS'] & fibermask.VARIABLETHRU) == 0)
+        skyfibers = np.where(oksky)[0]
         if exclude_sky_targetids is not None:
             log.info('Excluding default sky fibers using exclude_sky_targetids')
             bads = np.in1d(fibermap['TARGETID'][skyfibers], exclude_sky_targetids)

py/desispec/test/test_fiberbitmask.py

@@ -73,11 +73,10 @@ def test_ambiguous_maskbits(self):
         #- POORPOSITION is ok for flats, sky, and fluxcalib; but bad for stdstars
         self.check_mask('POORPOSITION', ok_steps=['flat', 'sky', 'fluxcalib'], bad_steps=['stdstar'])
 
-        #- NEARCHARGETRAP is informative; treated as ok for everyone including sky
+        #- NEARCHARGETRAP and VARIABLETHRU are informative for fiberbitmasking;
+        #- i.e. they don't trigger masking fibers
         #- TODO: it's actually bad for faint targets and sky for a single amp, but we structurally
         #- don't have a way to encode that in FIBERSTATUS (fiber not CCD or amp)
         self.check_mask('NEARCHARGETRAP', ok_steps=['flat', 'sky', 'stdstar', 'fluxcalib'], bad_steps=[])
+        self.check_mask('VARIABLETHRU', ok_steps=['flat', 'sky', 'stdstar', 'fluxcalib'], bad_steps=[])
 
-        #- VARIABLETHRU is ok for flats because otherwise we'd block the entire fiber,
-        #- and ok to at least attempt to flux calibrate it, but it shouldn't be used for sky or stdstars
-        self.check_mask('VARIABLETHRU', ok_steps=['flat', 'fluxcalib'], bad_steps=['sky', 'stdstar'])

py/desispec/test/test_flux_calibration.py

@@ -10,7 +10,7 @@
 import numpy as np
 #import scipy.sparse
 
-#from desispec.maskbits import specmask
+from desispec.maskbits import fibermask
 from desispec.frame import Frame
 from desispec.fluxcalibration import normalize_templates
 from desispec.fluxcalibration import FluxCalib
@@ -204,6 +204,7 @@ def test_isStdStar(self):
         fm['CMX_TARGET'] = np.zeros(10, dtype=int)
         fm['CMX_TARGET'][0:2] = cmx_mask.STD_FAINT
         fm['CMX_TARGET'][2:4] = cmx_mask.SV0_STD_FAINT
+        fm['FIBERSTATUS'] = 0
         self.assertEqual(main_cmx_or_sv(fm)[2], 'cmx')
         self.assertEqual(np.count_nonzero(isStdStar(fm)), 4)
 
@@ -213,6 +214,7 @@ def test_isStdStar(self):
         fm['SV1_MWS_TARGET'] = np.zeros(10, dtype=int)
         fm['SV1_DESI_TARGET'][0:2] = sv1_desi_mask.STD_FAINT
         fm['SV1_MWS_TARGET'][2:4] = sv1_mws_mask.GAIA_STD_FAINT
+        fm['FIBERSTATUS'] = 0
         self.assertEqual(main_cmx_or_sv(fm)[2], 'sv1')
         self.assertEqual(np.count_nonzero(isStdStar(fm)), 4)
 
@@ -224,11 +226,21 @@ def test_isStdStar(self):
         fm['DESI_TARGET'][2:4] = desi_mask.STD_BRIGHT
         fm['DESI_TARGET'][4:6] |= desi_mask.MWS_ANY
         fm['MWS_TARGET'][4:6] = sv1_mws_mask.GAIA_STD_FAINT
+        fm['FIBERSTATUS'] = 0
         self.assertEqual(main_cmx_or_sv(fm)[2], 'main')
         self.assertEqual(np.count_nonzero(isStdStar(fm)), 6)
         self.assertEqual(np.count_nonzero(isStdStar(fm, bright=False)), 4)
         self.assertEqual(np.count_nonzero(isStdStar(fm, bright=True)), 2)
 
+        #- VARIABLETHRU should be excluded
+        fm['FIBERSTATUS'] = 0
+        ii = isStdStar(fm)
+        fm['FIBERSTATUS'][0] = fibermask.VARIABLETHRU
+        jj = isStdStar(fm)
+        self.assertTrue(ii[0])
+        self.assertFalse(jj[0])
+
+
 
     def test_main(self):
         pass

py/desispec/test/test_sky.py

@@ -5,9 +5,11 @@
 import unittest
 
 import numpy as np
-from desispec.sky import compute_sky, subtract_sky, SkyModel
+from astropy.table import Table
+from desispec.sky import compute_sky, subtract_sky, SkyModel, get_sky_fibers
 from desispec.resolution import Resolution
 from desispec.frame import Frame
+from desispec.maskbits import fibermask
 import desispec.io
 
 import desispec.scripts.sky as skyscript
@@ -155,6 +157,29 @@ def test_sky_slice(self):
         self.assertEqual(sky2.throughput_corrections, None)
         self.assertEqual(sky2.nrej, sky1.nrej)
 
+    def test_get_sky_fibers(self):
+        """Test desispec.sky.get_sky_fibers"""
+        fibermap = Table()
+        fibermap['TARGETID'] = np.arange(5)
+        fibermap['OBJTYPE'] = ['TGT', 'SKY', 'TGT', 'SKY', 'TGT']
+        fibermap['FIBERSTATUS'] = 0
+
+        #- Test default behavior
+        skyfibers = list(get_sky_fibers(fibermap))
+        self.assertEqual(skyfibers, [1,3])
+
+        #- Test overrides
+        skyfibers = list(get_sky_fibers(fibermap, override_sky_targetids=[0,1]))
+        self.assertEqual(skyfibers, [0,1])
+
+        skyfibers = list(get_sky_fibers(fibermap, exclude_sky_targetids=[0,1]))
+        self.assertEqual(skyfibers, [3,])
+
+        #- SKY fibers with VARIABLETHRU should be excluded
+        fibermap['FIBERSTATUS'][3] = fibermask.VARIABLETHRU
+        skyfibers = list(get_sky_fibers(fibermap))
+        self.assertEqual(skyfibers, [1,])
+
     def test_main(self):
         pass