diff --git a/py/desisurvey/NTS.py b/py/desisurvey/NTS.py index 22fc7b8..75aad39 100644 --- a/py/desisurvey/NTS.py +++ b/py/desisurvey/NTS.py @@ -368,6 +368,16 @@ def next_tile(self, conditions=None, exposure=None, constraints=None, s2n = 50.0 * texp_remaining/texp_tot exptime = texp_remaining maxtime = self.ETC.MAX_EXPTIME + maxtimecond = getattr(self.config, 'maximum_time_in_conditions', + None) + if maxtimecond is not None: + maxtimecond = getattr(maxtimecond, sched_program, None) + if maxtimecond is not None: + maxtimecond = maxtimecond().to(u.day).value + if maxtimecond is None: + maxtimecond = np.inf + maxtime = np.min([maxtime, maxtimecond]) + days_to_seconds = 60*60*24 fivemin = 5/60/24 # 5 minutes... pretty arbitrary. if ((mjd <= self.scheduler.night_ephem['dusk']-fivemin) or diff --git a/py/desisurvey/optimize.py b/py/desisurvey/optimize.py index 65ca52d..49f35c5 100644 --- a/py/desisurvey/optimize.py +++ b/py/desisurvey/optimize.py @@ -134,6 +134,10 @@ def __init__(self, condition, lst_edges, lst_hist, subset=None, start=None, stop texp_nom = u.Quantity([ getattr(config.nominal_exposure_time, program)() for program in tiles.tileprogram[tile_sel]]) + moon_up_factor = getattr(config, 'moon_up_factor', None) + if moon_up_factor is not None: + moon_up_factor = getattr(moon_up_factor, condition)() + texp_nom *= moon_up_factor if completed is not None: completed = astropy.table.Table.read(completed) nobtained = np.zeros(tiles.ntiles, dtype='f4') @@ -146,6 +150,7 @@ def __init__(self, condition, lst_edges, lst_hist, subset=None, start=None, stop completed['NNIGHT_NEEDED_'+condition][mask]) boost = np.clip(nneeded-nobtained, 0, np.inf) texp_nom *= boost[tile_sel] + self.dlst_nom = 360 * texp_nom.to(u.day).value / 0.99726956583 # Save arrays for the tiles to plan. @@ -193,6 +198,10 @@ def __init__(self, condition, lst_edges, lst_hist, subset=None, start=None, stop self.use_plan(save_history=False) self.min_total_time = self.plan_hist.sum() + self.plan_tiles = self.get_plan(np.zeros(self.ntiles)) + self.use_plan() + self.plan_hist_ha0 = self.plan_hist.copy() + # Initialize HA assignments for each tile. if init == 'zero': self.ha = np.zeros(self.ntiles) @@ -220,6 +229,8 @@ def __init__(self, condition, lst_edges, lst_hist, subset=None, start=None, stop lst_cdf = np.zeros_like(edges) lst_cdf[1:] = np.cumsum(hist) lst_cdf /= lst_cdf[-1] + lst_cen = (edges[1:]+edges[:-1])/2 + lst_cdf = (lst_cdf[1:]+lst_cdf[:-1])/2 # Calculate the CDF of planned LST usage relative to the same # central LST, assuming HA=0. Instead of spreading each exposure # over multiple LST bins, add its entire HA=0 exposure time at @@ -229,9 +240,10 @@ def __init__(self, condition, lst_edges, lst_hist, subset=None, start=None, stop sort_idx = np.argsort(tile_ra) tile_cdf = np.cumsum(exptime[sort_idx]) tile_cdf /= tile_cdf[-1] + tile_cdf = (np.concatenate([[0], tile_cdf[:-1]])+tile_cdf)/2 # Use linear interpolation to find an LST for each tile that # matches the plan CDF to the available LST CDF. - new_lst = np.interp(tile_cdf, lst_cdf, edges) + new_lst = np.interp(tile_cdf, lst_cdf, lst_cen) # Calculate each tile's HA as the difference between its HA=0 # LST and its new LST after CDF matching. ha = np.empty(self.ntiles) diff --git a/py/desisurvey/scripts/surveyinit.py b/py/desisurvey/scripts/surveyinit.py index 8983bc8..8e9a258 100644 --- a/py/desisurvey/scripts/surveyinit.py +++ b/py/desisurvey/scripts/surveyinit.py @@ -98,6 +98,9 @@ def parse(options=None): parser.add_argument( '--completed', default=None, help='filename with information on already completed tiles') + parser.add_argument( + '--include-weather', default=True, type=bool, + help='Use past weather to discount available LST when planning.') if options is None: args = parser.parse_args() @@ -137,6 +140,8 @@ def calculate_initial_plan(args): for year in years: fractions.append( desimodel.weather.dome_closed_fractions(first, last, replay='Y{}'.format(year))) + from scipy.ndimage import gaussian_filter + fractions = gaussian_filter(fractions, 7, mode='wrap') weather = 1 - np.mean(fractions, axis=0) # Save the weather fractions as the primary HDU. hdr['FIRST'] = first.isoformat() @@ -152,8 +157,14 @@ def calculate_initial_plan(args): # Calculate the distribution of available LST in each program # during the nominal survey [start, stop). ilo, ihi = (start - first).days, (stop - first).days + print('weather', args.include_weather) + if args.include_weather: + tweather = weather[ilo:ihi] + else: + tweather = None lst_hist, lst_bins = ephem.get_available_lst( - nbins=args.nbins, weather=weather[ilo:ihi], include_twilight=args.include_twilight) + nbins=args.nbins, weather=tweather, + include_twilight=args.include_twilight) # Initialize the output results table. conditions = ['DARK', 'GRAY', 'BRIGHT'] @@ -227,6 +238,7 @@ def calculate_initial_plan(args): .format(condition, plan_sum, avail_sum, 1e2 * margin)) # Save planned LST usage. table['PLAN'] = opt.plan_hist + table['HA0'] = opt.plan_hist_ha0 hdus.append(fits.BinTableHDU(table, name=condition)) # Calculate exposure times in (solar) seconds.