compute dynamic saturation limit

measure_raw.py

@@ -261,7 +261,7 @@ def run(self, ps=None, focus=False, momentsize=5,
         camera = primhdr.get('INSTRUME','').strip().lower()
         # -> "decam" / "mosaic3"
         meas = dict(band=band, airmass=airmass, 
-                    skybright=skybr, primhdr=primhdr,
+                    skybright=skybr, pixscale=pixsc, primhdr=primhdr,
                     hdr=hdr, wcs=wcs, ra_ccd=ra_ccd, dec_ccd=dec_ccd,
                     extension=ext, camera=camera, 
                     ndetected=ndetected)

mosbot.py

@@ -27,9 +27,8 @@
 
 from mosaicstrategy import (
     ExposureFactor, getParserAndGlobals, setupGlobals,
-    GetAirmass, StartAndEndTimes, s_to_days, readTilesTable, GetNightlyStrategy,
-    WriteJSON)
-from measure_raw import measure_raw, get_default_extension
+    GetAirmass )
+from measure_raw import measure_raw, get_default_extension, get_nominal_cal
 
 from jnox import *
 from copilot import get_tile_from_name
@@ -484,14 +483,26 @@ def found_new_image(self, fn):
             exptime *= 1.1
             exptime = int(np.ceil(exptime))
             print('Exposure time with safety factor:', exptime)
-    
+
             exptime = np.clip(exptime, self.gvs.floor_exptimes[nextband],
                               self.gvs.ceil_exptimes[nextband])
             print('Clipped exptime', exptime)
-            if nextband == 'z' and exptime > self.gvs.t_sat_max:
-                exptime = self.gvs.t_sat_max
-                print('Reduced exposure time to avoid z-band saturation:',
-                      exptime)
+            if nextband == 'z':
+                # Compute cap on exposure time to avoid saturation /
+                # loss of dynamic range.
+                # Convert sky brightness in mag/arcsec^2 into counts
+
+                cal = get_nominal_cal(M['camera'], nextband)
+                nextzp0 = cal[0]
+
+                skyflux = 10. ** ((nextsky - nextzp0) / -2.5)
+                skyflux *= M['pixscale']**2
+                print('Predicted sky flux per pixel per second: %.1f' %skyflux)
+                t_sat = np.floor(30000. / skyflux)
+                if exptime > t_sat:
+                    exptime = t_sat
+                    print('Reduced exposure time to avoid z-band saturation:',
+                          exptime)
             exptime = int(exptime)
 
             print('Changing exptime from', jplan['expTime'], 'to', exptime)

test_mosbot.py

@@ -1,3 +1,4 @@
+from __future__ import print_function
 import unittest
 #from django.test import TestCase
 import os
@@ -14,7 +15,34 @@ def test_new_file(self):
         tempdir = tempfile.mkdtemp()
         args = ['--script', os.path.join(tempdir, 'tonight.sh')]
         ## FIXME -- will need to update these to have current dates; also put in git
-        args += ['pass1.json', 'pass2.json', 'pass3.json']
+        fn1 = os.path.join(self.testdatadir, 'pass1.json')
+        fn2 = os.path.join(self.testdatadir, 'pass2.json')
+        fn3 = os.path.join(self.testdatadir, 'pass3.json')
+
+        tmpfn1 = fn1 + '.tmp'
+        tmpfn2 = fn2 + '.tmp'
+        tmpfn3 = fn3 + '.tmp'
+
+        for fn,tmpfn in ((fn1,tmpfn1),(fn2,tmpfn2),(fn3,tmpfn3)):
+            import json
+            import ephem
+
+            J = json.loads(open(fn, 'r').read())
+            t0 = ephem.Date(str(J[0]['approx_datetime']))
+            print('First exposure:', t0)
+            now = ephem.now()
+            print('Now:', now)
+            for j in J:
+                tnew = now + ephem.Date(str(j['approx_datetime'])) - t0
+                tnew = str(ephem.Date(tnew))
+                j['approx_datetime'] = tnew
+                print('Updated datetime to', tnew)
+
+            f = open(tmpfn, 'w')
+            json.dump(J, f, sort_keys=True, indent=4, separators=(',', ': '))
+            f.close()
+
+        args += [tmpfn1, tmpfn2, tmpfn3]
 
         mosbot = main(cmdlineargs=args, get_mosbot=True)