|
| 1 | +#!/usr/bin/env python |
| 2 | + |
| 3 | +import os |
| 4 | +import fitsio |
| 5 | +import numpy as np |
| 6 | +from astropy.io import fits |
| 7 | +from astropy.table import Table, vstack |
| 8 | +from desiutil.log import get_logger |
| 9 | +from desiutil.redirect import stdouterr_redirected |
| 10 | +from desisurveyops.fba_tertiary_design_io import ( |
| 11 | + assert_environ_settings, |
| 12 | + get_fn, |
| 13 | + read_yaml, |
| 14 | + assert_tertiary_settings, |
| 15 | + get_tile_centers_grid, |
| 16 | + create_tiles_table, |
| 17 | + creates_priority_table, |
| 18 | + match_coord, |
| 19 | + print_samples_overlap, |
| 20 | + subsample_targets_avail, |
| 21 | + finalize_target_table, |
| 22 | + assert_files, |
| 23 | + create_targets_assign, |
| 24 | + plot_targets_assign, |
| 25 | +) |
| 26 | +from desitarget.targetmask import desi_mask, bgs_mask |
| 27 | +from argparse import ArgumentParser |
| 28 | + |
| 29 | +# AR message to ensure that some key settings in the environment are correctly set |
| 30 | +# AR put it at the very top, so that it appears first, and is not redirected |
| 31 | +# AR (=burried) in the log file |
| 32 | +assert_environ_settings() |
| 33 | + |
| 34 | +log = get_logger() |
| 35 | + |
| 36 | + |
| 37 | +def parse(): |
| 38 | + parser = ArgumentParser() |
| 39 | + parser.add_argument( |
| 40 | + "--yamlfn", |
| 41 | + help="path to the tertiary-config-PROGNUMPAD.yaml file", |
| 42 | + type=str, |
| 43 | + default=None, |
| 44 | + ) |
| 45 | + parser.add_argument( |
| 46 | + "--steps", |
| 47 | + help="comma-separated list of steps to execute (default='tiles,priorities,targets,run,diagnosis')", |
| 48 | + type=str, |
| 49 | + default="tiles,priorities,targets,run,diagnosis", |
| 50 | + ) |
| 51 | + parser.add_argument( |
| 52 | + "--dry_run", |
| 53 | + action="store_true", |
| 54 | + help="for the 'run' step, only print the commands on the prompt", |
| 55 | + ) |
| 56 | + parser.add_argument( |
| 57 | + "--log-stdout", |
| 58 | + "--log_stdout", |
| 59 | + action="store_true", |
| 60 | + help="log to stdout instead of redirecting to a file", |
| 61 | + ) |
| 62 | + args = parser.parse_args() |
| 63 | + for kwargs in args._get_kwargs(): |
| 64 | + log.info("{} = {}".format(kwargs[0], kwargs[1])) |
| 65 | + return args |
| 66 | + |
| 67 | + |
| 68 | +# AR purposefully offset the "usual" cosmos dec value by +0.3 deg |
| 69 | +# AR to avoid shallow ibis/m464 region |
| 70 | +def create_tiles(tileid_start, ntile, obsconds, outfn): |
| 71 | + tileids = np.arange(tileid_start, tileid_start + ntile, dtype=int) |
| 72 | + |
| 73 | + field_ra, field_dec = 150.1, 2.482 # AR 2.482 = 2.182 + 0.3 |
| 74 | + |
| 75 | + tileras, tiledecs = get_tile_centers_grid(field_ra, field_dec, npt=ntile) |
| 76 | + assert tileras.size == ntile |
| 77 | + assert tiledecs.size == ntile |
| 78 | + d = create_tiles_table(tileids, tileras, tiledecs, obsconds) |
| 79 | + d.write(outfn) |
| 80 | + |
| 81 | + |
| 82 | +# AR usual priority scheme (+1 when observed) |
| 83 | +def create_priorities(yamlfn, outfn): |
| 84 | + d = creates_priority_table(yamlfn) |
| 85 | + # AR print |
| 86 | + d.pprint_all() |
| 87 | + for sample in np.unique(d["TERTIARY_TARGET"]): |
| 88 | + sel = d["TERTIARY_TARGET"] == sample |
| 89 | + log.info("priorites for {}: {}".format(sample, d["PRIORITY"][sel].tolist())) |
| 90 | + d.write(outfn) |
| 91 | + |
| 92 | + |
| 93 | +def create_targets(yamlfn, outfn): |
| 94 | + |
| 95 | + # matching radius [arcsec] |
| 96 | + search_radius = 1 |
| 97 | + |
| 98 | + # AR read config |
| 99 | + config = read_yaml(yamlfn) |
| 100 | + prognum = config["settings"]["prognum"] |
| 101 | + targdir = config["settings"]["targdir"] |
| 102 | + rundate = config["settings"]["rundate"] |
| 103 | + |
| 104 | + samples = np.array(list(config["samples"].keys())) |
| 105 | + prios = np.array([config["samples"][sample]["PRIORITY_INIT"] for sample in samples]) |
| 106 | + ngoals = np.array([config["samples"][sample]["NGOAL"] for sample in samples]) |
| 107 | + |
| 108 | + # sort by *decreasing* priority |
| 109 | + ii = prios.argsort()[::-1] |
| 110 | + samples, prios, ngoals = samples[ii], prios[ii], ngoals[ii] |
| 111 | + log.info("reading {}:".format(yamlfn)) |
| 112 | + for sample in samples: |
| 113 | + msg = "{}:\t{}".format( |
| 114 | + sample, |
| 115 | + "\t".join( |
| 116 | + [ |
| 117 | + "{}={}".format(key, config["samples"][sample][key]) |
| 118 | + for key in config["samples"][sample] |
| 119 | + ] |
| 120 | + ), |
| 121 | + ) |
| 122 | + log.info(msg) |
| 123 | + log.info("") |
| 124 | + |
| 125 | + # read cats |
| 126 | + ds = {} |
| 127 | + log.info("reading cats from {}:".format(targdir)) |
| 128 | + for sample in samples: |
| 129 | + |
| 130 | + infn = os.path.join(targdir, "inputcats", config["samples"][sample]["FN"]) |
| 131 | + try: |
| 132 | + d = Table(fitsio.read(infn, columns=["RA", "DEC"])) |
| 133 | + except ValueError: |
| 134 | + d = Table(fitsio.read(infn, columns=["ra", "dec"])) |
| 135 | + d["PMRA"], d["PMDEC"], d["REF_EPOCH"] = 0.0, 0.0, 2015.5 |
| 136 | + |
| 137 | + d["ORIG_ROW"] = np.arange(len(d), dtype=int) |
| 138 | + d["ORIG_FN"] = infn |
| 139 | + d["CHECKER"] = config["samples"][sample]["CHECKER"] |
| 140 | + d["TERTIARY_TARGET"] = sample |
| 141 | + |
| 142 | + ds[sample] = d |
| 143 | + log.info("{}:\t{} rows".format(sample, len(d))) |
| 144 | + log.info("") |
| 145 | + |
| 146 | + # merge |
| 147 | + # AR merge |
| 148 | + log.info("merging (use {} arcsec search radius)".format(search_radius)) |
| 149 | + for sample, prio, ngoal in zip(samples, prios, ngoals): |
| 150 | + |
| 151 | + d = Table() |
| 152 | + for key in [ |
| 153 | + "RA", |
| 154 | + "DEC", |
| 155 | + "PMRA", |
| 156 | + "PMDEC", |
| 157 | + "REF_EPOCH", |
| 158 | + "TERTIARY_TARGET", |
| 159 | + "CHECKER", |
| 160 | + ]: |
| 161 | + d[key] = ds[sample][key] |
| 162 | + d["PRIORITY_INIT"] = prio |
| 163 | + d["NGOAL"] = ngoal |
| 164 | + log.info("{}:\tPRIORITY_INIT={}\tread {} rows".format(sample, prio, len(d))) |
| 165 | + |
| 166 | + for sample2 in samples: |
| 167 | + if sample2 == sample: |
| 168 | + d[sample2] = True |
| 169 | + d["{}_ROW".format(sample2)] = ds[sample]["ORIG_ROW"] |
| 170 | + else: |
| 171 | + d[sample2] = False |
| 172 | + d["{}_ROW".format(sample2)] = np.full(len(d), -99, dtype=int) |
| 173 | + |
| 174 | + if sample == samples[0]: |
| 175 | + log.info( |
| 176 | + "{}:\tPRIORITY_INIT={}\tadd all {} rows".format(sample, prio, len(d)) |
| 177 | + ) |
| 178 | + merge_d = d.copy() |
| 179 | + |
| 180 | + else: |
| 181 | + ii_merge, ii, _, _, _ = match_coord( |
| 182 | + merge_d["RA"], |
| 183 | + merge_d["DEC"], |
| 184 | + d["RA"], |
| 185 | + d["DEC"], |
| 186 | + search_radius=search_radius, |
| 187 | + verbose=False, |
| 188 | + ) |
| 189 | + # add info for matched ones |
| 190 | + log.info( |
| 191 | + "{}:\tPRIORITY_INIT={}\tmatched {} rows with targets ingested so far".format( |
| 192 | + sample, prio, ii.size |
| 193 | + ) |
| 194 | + ) |
| 195 | + merge_d[sample][ii_merge] = True |
| 196 | + merge_d["{}_ROW".format(sample)][ii_merge] = d["{}_ROW".format(sample)][ii] |
| 197 | + # add rows for non-matched ones |
| 198 | + sel = ~np.in1d(np.arange(len(d)), ii) |
| 199 | + d = d[sel] |
| 200 | + log.info( |
| 201 | + "{}:\tPRIORITY_INIT={}\tadd {} not-matched rows".format( |
| 202 | + sample, prio, len(d) |
| 203 | + ) |
| 204 | + ) |
| 205 | + merge_d = vstack([merge_d, d]) |
| 206 | + log.info("") |
| 207 | + |
| 208 | + # rename d, to make downstream code clearer |
| 209 | + d = merge_d.copy() |
| 210 | + |
| 211 | + # samples overlap |
| 212 | + log.info("looking at samples overlap before pre-filtering") |
| 213 | + print_samples_overlap(d, samples) |
| 214 | + log.info("") |
| 215 | + |
| 216 | + # AR finalize |
| 217 | + d = finalize_target_table(d, yamlfn) |
| 218 | + d.meta["RANDSEED"] = read_yaml(yamlfn)["settings"]["np_rand_seed"] |
| 219 | + d.meta["YAMLFN"] = yamlfn |
| 220 | + d.meta["SAMPLES"] = ",".join(samples) |
| 221 | + d.meta["PRIOS"] = ",".join(prios.astype(str)) |
| 222 | + d.meta["FNS"] = ",".join([config["samples"][sample]["FN"] for sample in samples]) |
| 223 | + d.write(outfn) |
| 224 | + |
| 225 | + |
| 226 | +def main(): |
| 227 | + |
| 228 | + # AR read + assert settings |
| 229 | + mydict = read_yaml(args.yamlfn)["settings"] |
| 230 | + assert_tertiary_settings(mydict) |
| 231 | + prognum, targdir = mydict["prognum"], mydict["targdir"] |
| 232 | + |
| 233 | + # AR set random seed (for SUBPRIORITY reproducibility) |
| 234 | + np.random.seed(mydict["np_rand_seed"]) |
| 235 | + |
| 236 | + # AR tiles file |
| 237 | + if "tiles" in args.steps.split(","): |
| 238 | + tilesfn = get_fn(prognum, "tiles", targdir) |
| 239 | + log.info("run create_tiles() to generate {}".format(tilesfn)) |
| 240 | + create_tiles( |
| 241 | + mydict["tileid_start"], mydict["ntile"], mydict["obsconds"], tilesfn |
| 242 | + ) |
| 243 | + |
| 244 | + # AR priorities file |
| 245 | + if "priorities" in args.steps.split(","): |
| 246 | + priofn = get_fn(prognum, "priorities", targdir) |
| 247 | + log.info("run create_priorities() to generate {}".format(priofn)) |
| 248 | + create_priorities(args.yamlfn, priofn) |
| 249 | + |
| 250 | + # AR targets file |
| 251 | + if "targets" in args.steps.split(","): |
| 252 | + targfn = get_fn(prognum, "targets", targdir) |
| 253 | + log.info("run create_targets() to generate {}".format(targfn)) |
| 254 | + create_targets(args.yamlfn, targfn) |
| 255 | + |
| 256 | + # AR sanity checks + run |
| 257 | + if "run" in args.steps.split(","): |
| 258 | + assert_files(prognum, targdir) |
| 259 | + cmd = "desi_fba_tertiary_wrapper --prognum {} --targdir {} --rundate {} --std_dtver {}".format( |
| 260 | + prognum, targdir, mydict["rundate"], mydict["std_dtver"] |
| 261 | + ) |
| 262 | + # AR required because of the safeguard... sigh.. |
| 263 | + # cmd = "{} --custom_too_development".format(cmd) |
| 264 | + if args.dry_run: |
| 265 | + cmd = "{} --dry_run".format(cmd) |
| 266 | + log.info(cmd) |
| 267 | + os.system(cmd) |
| 268 | + |
| 269 | + # AR diagnosis |
| 270 | + if "diagnosis" in args.steps.split(","): |
| 271 | + create_targets_assign(prognum, targdir) |
| 272 | + plot_targets_assign(prognum, targdir) |
| 273 | + |
| 274 | + |
| 275 | +if __name__ == "__main__": |
| 276 | + |
| 277 | + args = parse() |
| 278 | + |
| 279 | + if args.log_stdout: |
| 280 | + main() |
| 281 | + else: |
| 282 | + _ = read_yaml(args.yamlfn)["settings"] |
| 283 | + logfn = get_fn(_["prognum"], "log", _["targdir"]) |
| 284 | + with stdouterr_redirected(to=logfn): |
| 285 | + main() |
0 commit comments