Target selection for the first MWS Stream Program (GD1)

.github/workflows/python-package.yml

@@ -21,6 +21,7 @@ jobs:
                 astropy-version: ['==5.0', '<6']  # fuji version
                 fitsio-version: ['==1.1.6']  # fuji version
                 numpy-version: ['<1.23']  # to keep asscalar, used by astropy
+                matplotlib-version: ['<3.6.3'] # later versions of matplotlib require later versions of numpy.
 
         env:
             DESIUTIL_VERSION: 3.2.5
@@ -43,6 +44,7 @@ jobs:
                 python -m pip install -r requirements.txt
                 python -m pip install -U 'numpy${{ matrix.numpy-version }}'
                 python -m pip install -U 'astropy${{ matrix.astropy-version }}'
+                python -m pip install -U 'matplotlib${{ matrix.matplotlib-version }}'
                 python -m pip cache remove fitsio
                 python -m pip install --no-deps --force-reinstall --ignore-installed 'fitsio${{ matrix.fitsio-version }}'
                 svn export https://desi.lbl.gov/svn/code/desimodel/${DESIMODEL_DATA}/data
@@ -59,6 +61,7 @@ jobs:
                 python-version: [3.9]
                 fitsio-version: ['==1.1.6']  # fuji version
                 numpy-version: ['<1.23']  # to keep asscalar, used by astropy
+                matplotlib-version: ['<3.6.3'] # later versions of matplotlib require later versions of numpy.
         env:
             DESIUTIL_VERSION: 3.2.5
             DESIMODEL_DATA: branches/test-0.17
@@ -79,6 +82,7 @@ jobs:
                 python -m pip install git+https://github.com/desihub/desiutil.git@${DESIUTIL_VERSION}#egg=desiutil
                 python -m pip install -r requirements.txt
                 python -m pip install -U 'numpy${{ matrix.numpy-version }}'
+                python -m pip install -U 'matplotlib${{ matrix.matplotlib-version }}'
                 python -m pip cache remove fitsio
                 python -m pip install --no-deps --force-reinstall --ignore-installed 'fitsio${{ matrix.fitsio-version }}'
                 svn export https://desi.lbl.gov/svn/code/desimodel/${DESIMODEL_DATA}/data

bin/add_to_mtl_ledgers

@@ -0,0 +1,47 @@
+#!/usr/bin/env python
+
+from desitarget.mtl import add_to_ledgers, get_mtl_dir, _get_mtl_nside
+
+# ADM this is the upper-limit for memory for mtl._get_mtl_nside()=32.
+nproc = 12
+# ADM this is the upper-limit for memory for mtl._get_mtl_nside()=16.
+# nproc = 8
+# ADM retrieve the $MTL_DIR environment variable.
+mtldir = get_mtl_dir()
+mtlnside = _get_mtl_nside()
+# ADM default obsconditions.
+obscon = "DARK"
+
+from argparse import ArgumentParser
+ap = ArgumentParser(description='Add new targets to MTL ledger files from a file of targets')
+ap.add_argument("targfile",
+                help="Full path to a file of targets")
+ap.add_argument("--dest",
+                help="Full path to the output directory to host the ledger. The \
+                filename and full directory structure are built on-the-fly from \
+                file headers in targdirname. [defaults to {})".format(mtldir),
+                default=mtldir)
+ap.add_argument("--obscon",
+                help="String matching ONE obscondition in the bitmask yaml file \
+                (e.g. 'BRIGHT'). Controls priorities when merging targets and   \
+                where the ledger is written. [defaults to {})".format(obscon),
+                default=obscon)
+ap.add_argument("--numproc", type=int,
+                help='number of concurrent processes to use [defaults to {}]'.
+                format(nproc),
+                default=nproc)
+ap.add_argument('--healpixels',
+                help="HEALPixels corresponding to `nside` (e.g. '6,9,57'). Only \
+                write the ledger for targets in these pixels, at the default    \
+                nside, which is [{}]".format(mtlnside),
+                default=None)
+
+ns = ap.parse_args()
+
+# ADM parse the list of HEALPixels in which to run.
+pixlist = ns.healpixels
+if pixlist is not None:
+    pixlist = [int(pix) for pix in pixlist.split(',')]
+
+add_to_ledgers(ns.targfile, mtldir=ns.dest, pixlist=pixlist, obscon=ns.obscon,
+               numproc=ns.numproc)

bin/run_target_qa

@@ -79,7 +79,7 @@ else:
     max_bin_area = 1.0
 
 if ns.prepare:
-    fn = os.path.join(os.getenv("CSCRATCH"), ns.dest)
+    fn = os.path.join(os.getenv("SCRATCH"), ns.dest)
     objs, _, _, hdr = read_data(ns.src, header=True, downsample=ns.downsample)
     fitsio.write(fn+'.tmp', objs, extname='TARGETS', header=hdr, clobber=True)
     os.rename(fn+'.tmp', fn)

bin/select_stream_targets

@@ -0,0 +1,77 @@
+#!/usr/bin/env python
+
+import os, sys
+import numpy as np
+import fitsio
+
+from desitarget.streams import io, cuts
+
+from time import time
+start = time()
+
+from desiutil.log import get_logger
+log = get_logger()
+
+from argparse import ArgumentParser
+ap = ArgumentParser(description=("Generates DESI target bits from Legacy Surveys"
+                                 " sweep files")
+                    )
+ap.add_argument("sweepdir",
+                help=("Root directory of LS sweeps for a given data release for "
+                      "ONE of EITHER north or south")
+                )
+ap.add_argument("dest",
+                help=("Output target selection directory (the file name is built"
+                      " on-the-fly from other inputs)")
+                )
+ap.add_argument('-s','--streamnames', default="GD1",
+                help=("Comma-separated names of streams to run (e.g. x,y,z). "
+                      "Default is to just run GD1")
+                )
+ap.add_argument("--donotaddnors", action="store_true",
+                help=("Both the south/north LS files are read by default. Pass "
+                      "this to read only the specific files in passed sweepdir")
+                )
+ap.add_argument("--donotreadperstream", action="store_true",
+                help=("Default is to read targets in a loop per-stream instead "
+                      "of in a loor per-sweeps file. Pass this to loop over "
+                      " sweep files instead of streams")
+                )
+ap.add_argument("--donotreadcache", action="store_true",
+                help=("Default is to read from previously cached data files "
+                      "where possible. Pass this to start from scratch (and "
+                      "overwrite caches). Only relevant for --readperstream")
+                )
+
+ns = ap.parse_args()
+
+# ADM build the list of command line arguments to
+# ADM write to the output header.
+hdr = fitsio.FITSHDR()
+nsdict = vars(ns)
+for k in nsdict:
+    hdr[k.upper()] = nsdict[k]
+# ADM also explicitly add the command used, just in case.
+hdr["CMDLINE"] = ' '.join(sys.argv)
+
+# ADM change the input comma-separated string of stream names to a list.
+sn = ns.streamnames
+if sn is not None:
+    streamnames = sn.split(',')
+
+# ADM flip the sense of some inputs.
+addnors = not(ns.donotaddnors)
+readperstream = not(ns.donotreadperstream)
+readcache = not(ns.donotreadcache)
+
+targets = cuts.select_targets(
+    ns.sweepdir, stream_names=streamnames, readperstream=readperstream,
+    addnors=addnors, readcache=readcache
+)
+
+# ADM note that we hardcode for bright conditions. This could change in
+# ADM the future if we add new streams or programs.
+ntargs, outfile = io.write_targets(ns.dest, targets, hdr,
+                                   streamnames=sn, obscon="BRIGHT")
+
+log.info(f"{ntargs} targets written to {outfile}...t={time()-start:.1f}s")

doc/changes.rst

@@ -13,7 +13,13 @@ desitarget Change Log
     * Addresses `issue #818`_.
 * Function to match RA/Dec positions to Main Survey targets [`PR #820`_].
 * Bump astropy from 5.0 to 5.3.3 (dependabot) [`PR #815`_].
+* Add the GD1 dark matter streams program [`PR #814`_].
+  * Substantial new code in the `streams` directory.
+  * Addresses `issue #812`_.
+  * General unit tests fixes, including one introduced in `PR #823`_
 
+.. _`issue #812`: https://github.com/desihub/desitarget/issues/812
+.. _`PR #814`: https://github.com/desihub/desitarget/pull/814
 .. _`PR #815`: https://github.com/desihub/desitarget/pull/815
 .. _`issue #818`: https://github.com/desihub/desitarget/issues/818
 .. _`PR #820`: https://github.com/desihub/desitarget/pull/820

py/desitarget/cuts.py

@@ -2392,7 +2392,7 @@ def set_target_bits(photsys_north, photsys_south, obs_rflux,
     south_cuts = [False, True]
     if resolvetargs:
         # ADM if only southern objects were sent this will be [True], if
-        # ADM only northern it will be [False], else it wil be both.
+        # ADM only northern it will be [False], else it will be both.
         south_cuts = list(set(np.atleast_1d(photsys_south)))
 
     # ADM default for target classes we WON'T process is all False.
@@ -2599,13 +2599,21 @@ def set_target_bits(photsys_north, photsys_south, obs_rflux,
 
     # ADM when resolving, strictly only set bits for targets that pass
     # ADM the northern (southern) cuts in northern (southern) imaging.
-    res_north = np.ones_like(photsys_north)
-    res_south = np.ones_like(photsys_south)
-    if resolvetargs:
-        # ADM this construction guards against mutability but is far
-        # ADM quicker than making a copy.
-        res_north[:] = photsys_north[:]
-        res_south[:] = photsys_south[:]
+    # ADM first guard against photsys quantities bring passed as scalars.
+    if isinstance(photsys_north, bool):
+        res_north = True
+        res_south = True
+        if resolvetargs:
+            res_north = photsys_north
+            res_south = photsys_south
+    else:
+        res_north = np.ones_like(photsys_north)
+        res_south = np.ones_like(photsys_south)
+        if resolvetargs:
+            # ADM this construction guards against mutability but is far
+            # ADM quicker than making a copy.
+            res_north[:] = photsys_north[:]
+            res_south[:] = photsys_south[:]
 
     # Construct the targetflag bits for DECaLS (i.e. South).
     desi_target = (lrg_south & res_south) * desi_mask.LRG_SOUTH
@@ -2621,7 +2629,7 @@ def set_target_bits(photsys_north, photsys_south, obs_rflux,
     desi_target |= (elg_lop_north & res_north) * desi_mask.ELG_LOP_NORTH
     desi_target |= (qso_north & res_north) * desi_mask.QSO_NORTH
 
-    # Construct the targetflag bits combining north and south.
+    # Construct the target flag bits combining north and south.
     desi_target |= lrg * desi_mask.LRG
     desi_target |= elg * desi_mask.ELG
     desi_target |= elg_vlo * desi_mask.ELG_VLO

py/desitarget/data/streams.yaml

@@ -0,0 +1,83 @@
+# ADM This yaml file stores isochrones and coordinate information for stellar streams.
+
+# ADM the gd1 stream.
+GD1:
+    # ADM stream coordinates.
+    RAPOL: 34.5987
+    DECPOL: 29.7331
+    RA_REF: 200
+    # ADM stream extent.
+    MIND: 80
+    MAXD: 100
+    # ADM color offset.
+    COLOFF: -0.01
+    # ADM magnitude offset.
+    MAGOFF: -0.1
+    # ADM g-band isochrone.
+    ISO_G: [
+        14.1916, 13.8004, 13.4031, 12.9995, 12.6317, 12.3286, 12.0868, 11.8791,
+        11.6837, 11.4824, 11.2802, 11.0934, 10.9082, 10.7581, 10.6532, 10.586,
+        10.5413, 10.5083, 10.4813, 10.4547, 10.426, 10.3941, 10.3567, 10.3109,
+        10.2414, 9.9928, 9.7144, 9.5754, 9.4487, 9.3196, 9.1814, 9.0267,
+        8.8643, 8.706, 8.5551, 8.4145, 8.2813, 8.15, 8.0148, 7.8739, 7.7934,
+        7.7147, 7.637, 7.5598, 7.4838, 7.4089, 7.3364, 7.2656, 7.1963, 7.1287,
+        7.0627, 6.9984, 6.9357, 6.8747, 6.8153, 6.7572, 6.7005, 6.645, 6.5907,
+        6.5374, 6.4851, 6.4338, 6.3837, 6.3347, 6.2869, 6.2399, 6.1939, 6.1487,
+        6.1045, 6.0611, 5.9994, 5.9382, 5.8771, 5.8163, 5.7558, 5.6955, 5.6354,
+        5.5754, 5.5156, 5.4559, 5.3965, 5.3375, 5.2791, 5.2213, 5.1642, 5.1079,
+        5.0524, 4.9977, 4.9436, 4.8902, 4.8372, 4.7848, 4.7329, 4.6817, 4.6309,
+        4.5808, 4.5312, 4.482, 4.4333, 4.385, 4.3699, 4.3544, 4.339, 4.3232,
+        4.3073, 4.2913, 4.2752, 4.259, 4.2426, 4.2262, 4.2096, 4.1929, 4.176,
+        4.159, 4.1419, 4.1246, 4.1073, 4.0899, 4.0723, 4.0547, 4.0369, 4.0192,
+        4.0013, 3.9833, 3.9652, 3.9471, 3.9289, 3.9106, 3.8923, 3.8739, 3.8555,
+        3.837, 3.8185, 3.7999, 3.7814, 3.7629, 3.7444, 3.726, 3.7077, 3.6895,
+        3.6715, 3.6537, 3.6362, 3.6189, 3.6018, 3.5851, 3.5688, 3.5529, 3.5373,
+        3.5228, 3.5145, 3.5072, 3.5006, 3.4933, 3.4841, 3.4738, 3.4632, 3.4515,
+        3.4384, 3.4234, 3.4063, 3.387, 3.3656, 3.3425, 3.3181, 3.2925, 3.2661,
+        3.2391, 3.2118, 3.184, 3.1526, 3.1202, 3.0872, 3.0534, 3.0192, 2.9844,
+        2.9492, 2.9137, 2.8777, 2.8416, 2.805, 2.7678, 2.7301, 2.6921, 2.6533,
+        2.6143, 2.5749, 2.5341, 2.4911, 2.4469, 2.4014, 2.3549, 2.3073, 2.2589,
+        2.2096, 2.1599, 2.1095, 2.0585, 2.0069, 1.9546, 1.9018, 1.8486, 1.7948,
+        1.7407, 1.6863, 1.6768, 1.6668, 1.6566, 1.6464, 1.6365, 1.6264, 1.6162,
+        1.6064, 1.5965, 1.5863, 1.474, 1.3721, 1.2678, 1.1626, 1.063, 0.9624,
+        0.8616, 0.7646, 0.6683, 0.5721, 0.4812, 0.3935, 0.3142, 0.2533, 0.2509,
+        0.2242, 0.1251, 0.0243, -0.0714, -0.165, -0.2543, -0.3415, -0.4262,
+        -0.5093, -0.5895, -0.6659, -0.7401, -0.8117, -0.8808, -0.9467, -1.0111,
+        -1.0726, -1.1319, -1.1889, -1.2441, -1.2969, -1.3477, -1.3951, -1.4404,
+        -1.4833, -1.5249, -1.5657, -1.6049, -1.6429, -1.6793, -1.7145, -1.7484,
+        -1.7809, -1.8122, -1.8428, -1.8719, -1.8998, -1.9234
+    ]
+    # ADM r-band isochrone.
+    ISO_R: [
+        12.5842, 12.2175, 11.8471, 11.4737, 11.1368, 10.8625, 10.6455, 10.4599,
+        10.2854, 10.1058, 9.9258, 9.7602, 9.5962, 9.4637, 9.3714, 9.3122,
+        9.273, 9.2441, 9.2205, 9.1974, 9.1726, 9.1447, 9.112, 9.0719, 9.0111,
+        8.7983, 8.5647, 8.4486, 8.3433, 8.2368, 8.1238, 7.9988, 7.8693, 7.7443,
+        7.6265, 7.5175, 7.4151, 7.3147, 7.2116, 7.1038, 7.0416, 6.9804, 6.9199,
+        6.86, 6.8008, 6.7426, 6.6853, 6.629, 6.5737, 6.5194, 6.4663, 6.4143,
+        6.3635, 6.3138, 6.2652, 6.2177, 6.1712, 6.1257, 6.081, 6.0372, 5.9942,
+        5.952, 5.9105, 5.8698, 5.8298, 5.7905, 5.7517, 5.7137, 5.6763, 5.6395,
+        5.587, 5.5345, 5.4821, 5.4298, 5.3774, 5.3252, 5.273, 5.2211, 5.1692,
+        5.1172, 5.0653, 5.0136, 4.9622, 4.9112, 4.8607, 4.8107, 4.761, 4.7117,
+        4.663, 4.6148, 4.5673, 4.5202, 4.4734, 4.4268, 4.3804, 4.3341, 4.2879,
+        4.2418, 4.196, 4.1503, 4.1361, 4.1216, 4.1071, 4.0924, 4.0777, 4.0628,
+        4.0478, 4.0327, 4.0174, 4.0018, 3.986, 3.97, 3.9538, 3.9374, 3.9208,
+        3.904, 3.8871, 3.87, 3.8528, 3.8354, 3.8178, 3.8001, 3.782, 3.7638,
+        3.7454, 3.7268, 3.708, 3.689, 3.6697, 3.6502, 3.6306, 3.6108, 3.5909,
+        3.5707, 3.5504, 3.5299, 3.5093, 3.4885, 3.4675, 3.4464, 3.425, 3.4035,
+        3.3819, 3.3601, 3.3382, 3.3163, 3.2942, 3.2721, 3.2501, 3.2281, 3.2084,
+        3.1879, 3.1675, 3.1473, 3.1269, 3.1062, 3.0853, 3.0637, 3.0411, 3.0169,
+        2.9913, 2.9644, 2.9363, 2.9075, 2.8779, 2.8478, 2.8174, 2.7866, 2.7555,
+        2.7242, 2.6892, 2.6536, 2.6175, 2.5809, 2.5442, 2.5071, 2.4696, 2.432,
+        2.3942, 2.3563, 2.318, 2.2792, 2.2401, 2.2004, 2.1602, 2.1197, 2.079,
+        2.0367, 1.9924, 1.9466, 1.8995, 1.8514, 1.8021, 1.752, 1.7009, 1.6492,
+        1.5967, 1.5434, 1.4894, 1.4348, 1.3796, 1.3238, 1.2675, 1.2106, 1.1533,
+        1.1429, 1.1325, 1.1219, 1.1113, 1.1007, 1.0901, 1.0794, 1.0687, 1.058,
+        1.0472, 0.9282, 0.8196, 0.7079, 0.5949, 0.487, 0.3786, 0.2698, 0.165,
+        0.0607, -0.0441, -0.1435, -0.2399, -0.3271, -0.394, -0.3956, -0.4237,
+        -0.5341, -0.6472, -0.7552, -0.8614, -0.9633, -1.0632, -1.1607, -1.2567,
+        -1.3499, -1.4398, -1.5283, -1.6146, -1.6986, -1.7795, -1.8595, -1.9366,
+        -2.0117, -2.0847, -2.1561, -2.2251, -2.2923, -2.357, -2.42, -2.4803,
+        -2.5395, -2.5976, -2.6541, -2.7092, -2.7626, -2.8146, -2.8653, -2.9145,
+        -2.9622, -3.0091, -3.0542, -3.0979, -3.1371
+    ]

py/desitarget/data/targetmask.yaml

@@ -174,6 +174,9 @@ scnd_mask:
     - [DWF_BRIGHT_LO,          52, "See $SCND_DIR/docs/DWF_BRIGHT_LO.txt",        {obsconditions: BRIGHT,      filename: 'DWF_BRIGHT_LO',       flavor: 'SPARE', updatemws: True,  downsample: 1}]
     - [DWF_DARK_HI,            53, "See $SCND_DIR/docs/DWF_DARK_HI.txt",          {obsconditions: DARK,        filename: 'DWF_DARK_HI',         flavor: 'SPARE', updatemws: True,  downsample: 1}]
     - [DWF_DARK_LO,            54, "See $SCND_DIR/docs/DWF_DARK_LO.txt",          {obsconditions: DARK,        filename: 'DWF_DARK_LO',         flavor: 'SPARE', updatemws: True,  downsample: 1}]
+    - [GD1_BRIGHT_PM,          55, "Bright targets with Gaia in GD1 stream",      {obsconditions: BRIGHT,      filename: None,                  flavor: 'SPARE', updatemws: True,  downsample: 1}]
+    - [GD1_FAINT_NO_PM,        56, "Faint targets without Gaia in GD1 stream",    {obsconditions: BRIGHT,      filename: None,                  flavor: 'SPARE', updatemws: True,  downsample: 1}]
+    - [GD1_FILLER,             57, "Filler targets in GD1 stream",                {obsconditions: BRIGHT,      filename: None,                  flavor: 'SPARE', updatemws: True,  downsample: 1}]
 
 
 # ADM reserve 59-62 in scnd_mask for Targets of Opportunity in both SV and the Main Survey.
@@ -414,6 +417,9 @@ priorities:
         BRIGHT_TOO_HIP:         {UNOBS: 9500, MORE_ZGOOD: 2,    MORE_ZWARN: 2,    DONE: 2, OBS: 1, DONOTOBSERVE: 0, MORE_MIDZQSO: 0}
         DARK_TOO_LOP:           SAME_AS_BRIGHT_TOO_LOP
         DARK_TOO_HIP:           SAME_AS_BRIGHT_TOO_HIP
+        GD1_BRIGHT_PM:          {UNOBS: 1521, MORE_ZGOOD: 1520, MORE_ZWARN: 1520, DONE: 2, OBS: 1, DONOTOBSERVE: 0, MORE_MIDZQSO: 0}
+        GD1_FAINT_NO_PM:        {UNOBS: 1510, MORE_ZGOOD: 1511, MORE_ZWARN: 1511, DONE: 2, OBS: 1, DONOTOBSERVE: 0, MORE_MIDZQSO: 0}
+        GD1_FILLER:             {UNOBS: 49,   MORE_ZGOOD: 50,   MORE_ZWARN: 50,   DONE: 2, OBS: 1, DONOTOBSERVE: 0, MORE_MIDZQSO: 0}
 
 # ADM INITIAL number of observations (NUMOBS_INIT) for each target bit
 # ADM SAME_AS_XXX means to use the NUMOBS_INIT for bitname XXX
@@ -565,3 +571,6 @@ numobs:
         DWF_BRIGHT_LO:          10
         DWF_DARK_HI:            10
         DWF_DARK_LO:            10
+        GD1_BRIGHT_PM:          7
+        GD1_FAINT_NO_PM:        7
+        GD1_FILLER:             7