Merge pull request #1123 from desihub/darks

darkmodel bookkeeping; desi_compute_psf --specmin bugfix

bin/desi_compute_bias

@@ -8,13 +8,20 @@ parser = argparse.ArgumentParser(formatter_class=argparse.ArgumentDefaultsHelpFo
 description="Compute a master bias from a set of raw data bias images",
 epilog='''This is simply a median of the input raw images.'''
 )
-parser.add_argument('-i','--image', type = str, default = None, required = True, nargs="*",
+parser.add_argument('-i','--image', type = str, default = None, nargs="*",
                     help = 'path of image fits files')
-parser.add_argument('-o','--outfile', type = str, default = None, required = True,
+parser.add_argument('-o','--outfile', type=str, default = None, required = True,
                     help = 'output median image filename')
 parser.add_argument('--camera',type = str, required = True,
                     help = 'camera name BX,RX,ZX with X from 0 to 9')
+parser.add_argument('--explistfile', type=str, default=None, required=False,
+                    help = 'Read NIGHT EXPID from file (one per line)')
 
 args        = parser.parse_args()
 
-compute_bias_file(args.image, args.outfile, args.camera)
+if (args.image is not None) and (args.explistfile is not None):
+    print('ERROR: specify --image or --explistfile but not both')
+    sys.exit(1)
+
+compute_bias_file(args.image, args.outfile, args.camera,
+                  explistfile=args.explistfile)

bin/desi_compute_dark_nonlinear

@@ -18,10 +18,17 @@ parser = argparse.ArgumentParser(
     model(x,y,t) = bias(x,y) + dark(x,y)*t + nonlinear(y,t)
 
     i.e. the non-linear term is only a function of row (y)
+
+    Data can be grouped by calendar --days (rollover at midnight) or
+    --nights (rollover at noon)
 ''')
 
-parser.add_argument('--days', type=int, required=True, nargs="*",
+parser.add_argument('--days', type=int, nargs="*",
                     help='YEARMMDD days to use for ZEROs and DARKs')
+parser.add_argument('--nights', type=int, nargs="*",
+                    help='YEARMMDD nights to use for ZEROs and DARKs')
+parser.add_argument('--first-expid', type=int, required=False,
+                    help='First EXPID to include')
 parser.add_argument('--darkfile', type=str, required=True,
                     help='output dark model file')
 parser.add_argument('--biasfile', type=str, required=True,
@@ -59,6 +66,10 @@ from desispec.io import findfile
 
 log = get_logger()
 
+if args.days is None and args.nights is None:
+    log.critical('Must specify --days or --nights')
+    sys.exit(1)
+
 #- tempdir caches files that could be re-used when rerunning for debugging
 #- i.e. it can be cleaned up when done, but isn't completely transient
 if args.tempdir is None:
@@ -73,13 +84,21 @@ if not os.path.isdir(tempdir):
 
 #- Data taken on the morning of the first day might be associated with the
 #- previous NIGHT
-year = int(str(args.days[0])[0:4])
-month = int(str(args.days[0])[4:6])
-day = int(str(args.days[0])[6:8])
-
-t = datetime.datetime(year, month, day) - datetime.timedelta(days=1)
-nights = [int(t.strftime('%Y%m%d'))]
-nights.extend(args.days)
+#- Note: use the variable named "night" whether grouped by days or nights
+if args.days:
+    log.info('Grouping ZEROs with DARKs by calendar DAY')
+    daynight = 'DAY'
+    year = int(str(args.days[0])[0:4])
+    month = int(str(args.days[0])[4:6])
+    day = int(str(args.days[0])[6:8])
+
+    t = datetime.datetime(year, month, day) - datetime.timedelta(days=1)
+    nights = [int(t.strftime('%Y%m%d'))]
+    nights.extend(args.days)
+else:
+    log.info('Grouping ZEROs with DARKs by observing NIGHT')
+    daynight = 'NIGHT'
+    nights = args.nights
 
 #- Get table of what exposures were taken on those days
 speclog_file = os.path.join(tempdir, 'speclog.csv')
@@ -91,20 +110,29 @@ else:
     speclog = desispec.io.util.get_speclog(nights)
 
     #- Add "DAY" column = rolls over at midnight instead of MST noon
-    t = Time(speclog['MJD']-7/24, format='mjd')
-    speclog['DAY'] = t.strftime('%Y%m%d').astype(int)
+    if daynight == 'DAY':
+        t = Time(speclog['MJD']-7/24, format='mjd')
+        speclog['DAY'] = t.strftime('%Y%m%d').astype(int)
 
-    #- Trim to just the requested days
+    #- Trim to just the requested days/nights
     keep = np.zeros(len(speclog), dtype=bool)
-    for day in args.days:
-        keep |= speclog['DAY'] == day
+    for night in nights:
+        keep |= speclog[daynight] == night
 
     speclog = speclog[keep]
     tmpfile = speclog_file + '.tmp-' + str(os.getpid())
     speclog.write(tmpfile, format='ascii.csv')
     os.rename(tmpfile, speclog_file)
     log.info(f'Wrote speclog to {speclog_file}')
 
+#- filter expids if requested
+if args.first_expid is not None:
+    keep = speclog['EXPID'] >= args.first_expid
+    speclog = speclog[keep]
+    if len(speclog) == 0:
+        log.critical('No exposures! exiting...')
+        sys.exit(1)
+
 #- group EXPTIMEs by integer
 speclog['EXPTIME_INT'] = speclog['EXPTIME'].astype(int)
 
@@ -124,35 +152,35 @@ speclog = speclog[keep]
 #- Print some summary stats before continuing
 isZero = speclog['OBSTYPE'] == 'ZERO'
 isDark = speclog['OBSTYPE'] == 'DARK'
-for day in args.days:
-    ii = speclog['DAY'] == day
+for night in nights:
+    ii = speclog[daynight] == night
     nzeros = np.count_nonzero(ii & isZero)
     ndarks = np.count_nonzero(ii & isDark)
     darktimes = sorted(set(speclog['EXPTIME_INT'][ii & isDark]))
-    log.info(f'Day {day} has {nzeros} ZEROs and {ndarks} DARKs with exptimes {darktimes}')
+    log.info(f'{night} has {nzeros} ZEROs and {ndarks} DARKs with exptimes {darktimes}')
 
 #- Combine the ZEROs into per-day bias files
 all_zerofiles = list()
-for day in args.days:
+for night in nights:
     zerofiles = list()
-    ii = isZero & (speclog['DAY'] == day)
+    ii = isZero & (speclog[daynight] == night)
     nzeros = np.count_nonzero(ii)
     nzeros_good = nzeros - args.nskip_zeros
     if nzeros_good < 5:
-        log.critical(f'{nzeros} ZEROS on {day} is insufficient when skipping {args.nskip_zeros}')
+        log.critical(f'{nzeros} ZEROS on {night} is insufficient when skipping {args.nskip_zeros}')
         sys.exit(1)
 
     elif nzeros_good < 20:
-        log.warning(f'Only {nzeros_good} good ZEROs on day {day}')
+        log.warning(f'Only {nzeros_good} good ZEROs on {night}')
     else:
-        log.info(f'Using {nzeros_good} ZEROs on day {day}')
+        log.info(f'Using {nzeros_good} ZEROs on {night}')
 
     for row in speclog[ii][args.nskip_zeros:]:
         rawfile = findfile('raw', row['NIGHT'], row['EXPID'])
         zerofiles.append(rawfile)
         all_zerofiles.append(rawfile)
 
-    biasfile = f'{tempdir}/bias-{day}-{args.camera}.fits'
+    biasfile = f'{tempdir}/bias-{night}-{args.camera}.fits'
     if os.path.exists(biasfile):
         log.info(f'{biasfile} already exists')
     else:
@@ -180,7 +208,7 @@ for exptime in darktimes:
     biasfiles = list()
     ii = (speclog['EXPTIME_INT'] == exptime)
     for row in speclog[isDark & ii]:
-        day, night, expid = row['DAY'], row['NIGHT'], row['EXPID']
+        day, night, expid = row[daynight], row['NIGHT'], row['EXPID']
         rawfiles.append(findfile('raw', night, expid, args.camera))
         biasfiles.append(f'{tempdir}/bias-{day}-{args.camera}.fits')
 

py/desispec/ccdcalib.py

@@ -187,17 +187,44 @@ def compute_dark_file(rawfiles, outfile, camera, bias=None, nocosmic=False,
     log.info(f"done")
 
 
-def compute_bias_file(rawfiles, outfile, camera):
+def compute_bias_file(rawfiles, outfile, camera, explistfile=None):
     """
     Compute a bias file from input ZERO rawfiles
 
     Args:
         rawfiles: list of input raw file names
         outfile (str): output filename
         camera (str): camera, e.g. b0, r1, z9
+
+    Options:
+        explistfile: filename with text list of NIGHT EXPID to use
+
+    Notes: explistfile is only used if rawfiles=None; it should have
+    one NIGHT EXPID entry per line.
     """
     log = get_logger()
 
+    if explistfile is not None:
+        if rawfiles is not None:
+            msg = "specify rawfiles or explistfile, but not both"
+            log.error(msg)
+            raise ValueError(msg)
+
+        rawfiles = list()
+        with open(explistfile, 'r') as fx:
+            for line in fx:
+                line = line.strip()
+                if line.startswith('#') or len(line)<2:
+                    continue
+                night, expid = map(int, line.split())
+                filename = io.findfile('raw', night, expid)
+                if not os.path.exists(filename):
+                    msg = f'Missing {filename}'
+                    log.critical(msg)
+                    raise RuntimeError(msg)
+
+                rawfiles.append(filename)
+
     log.info("read images ...")
     images=[]
     shape=None
@@ -363,10 +390,25 @@ def model_y1d(image, smooth=0):
 
     return model1d
 
-def make_dark_scripts(days, outdir, cameras=None,
-        linexptime=None, nskip_zeros=None, tempdir=None, nosubmit=False):
+def make_dark_scripts(outdir, days=None, nights=None, cameras=None,
+        linexptime=None, nskip_zeros=None, tempdir=None, nosubmit=False,
+        first_expid=None):
     """
-    TODO: document
+    Generate batch script to run desi_compute_dark_nonlinear
+
+    Args:
+        outdir (str): output directory
+        days or nights (list of int): days or nights to include
+
+    Options:
+        cameras (list of str): cameras to include, e.g. b0, r1, z9
+        linexptime (float): exptime after which dark current is linear
+        nskip_zeros (int): number of ZEROs at beginning of day/night to skip
+        tempdir (str): tempfile working directory
+        nosubmit (bool): generate scripts but don't submit them to batch queue
+        first_expid (int): ignore expids prior to this
+
+    Args/Options are passed to the desi_compute_dark_nonlinear script
     """
     log = get_logger()
 
@@ -389,7 +431,14 @@ def make_dark_scripts(days, outdir, cameras=None,
     today = datetime.datetime.now().strftime('%Y%m%d')
 
     #- Convert list of days to single string to use in command
-    days = ' '.join([str(tmp) for tmp in days])
+    if days is not None:
+        days = ' '.join([str(tmp) for tmp in days])
+    elif nights is not None:
+        nights = ' '.join([str(tmp) for tmp in nights])
+    else:
+        msg = 'Must specify days or nights'
+        log.critical(msg)
+        raise ValueError(msg)
 
     for camera in cameras:
         sp = 'sp' + camera[1]
@@ -400,15 +449,21 @@ def make_dark_scripts(days, outdir, cameras=None,
         darkfile = f'dark-{key}.fits.gz'
         biasfile = f'bias-{key}.fits.gz'
 
-        cmd = f"desi_compute_dark_nonlinear --days {days}"
+        cmd = f"desi_compute_dark_nonlinear"
         cmd += f" \\\n    --camera {camera}"
         cmd += f" \\\n    --tempdir {tempdir}"
         cmd += f" \\\n    --darkfile {darkfile}"
         cmd += f" \\\n    --biasfile {biasfile}"
+        if days is not None:
+            cmd += f" \\\n    --days {days}"
+        if nights is not None:
+            cmd += f" \\\n    --nights {nights}"
         if linexptime is not None:
             cmd += f" \\\n    --linexptime {linexptime}"
         if nskip_zeros is not None:
             cmd += f" \\\n    --nskip-zeros {nskip_zeros}"
+        if first_expid is not None:
+            cmd += f" \\\n    --first-expid {first_expid}"
 
         with open(batchfile, 'w') as fx:
             fx.write(f'''#!/bin/bash -l

py/desispec/scripts/specex.py

@@ -156,13 +156,12 @@ def main(args, comm=None):
 
 
     mynbundle = int(nbundle / nproc)
-    myfirstbundle = 0
     leftover = nbundle % nproc
     if rank < leftover:
         mynbundle += 1
-        myfirstbundle = rank * mynbundle
+        myfirstbundle = bundles[0] + rank * mynbundle
     else:
-        myfirstbundle = ((mynbundle + 1) * leftover) + \
+        myfirstbundle = bundles[0] + ((mynbundle + 1) * leftover) + \
             (mynbundle * (rank - leftover))
 
     if rank == 0: