diff --git a/pipeline/report.py b/pipeline/report.py new file mode 100644 index 0000000..2b1b025 --- /dev/null +++ b/pipeline/report.py @@ -0,0 +1,188 @@ +"""Stage 6: write the session's own METHODS.md. + +Every session directory gets an account of what was done to it and what came +out, written for someone who was not there and has no access to the chat that +produced it. That is the whole point of generating it rather than writing it by +hand: a hand-written note is written once and then rots, while this is rebuilt +from the actual results every time the pipeline runs. + +It records what was measured, not what was intended - including the things that +did not work, because a session where the plate solve failed is more useful +with that stated than with it omitted. +""" +import os +from datetime import datetime, timezone + +import numpy as np + +import layout + + +def _fmt_exposure(seconds): + return f"{seconds / 60:.0f} min" if seconds >= 60 else f"{seconds:.0f} s" + + +def write(session, results, verbose=True): + """Compose METHODS.md from whatever the stages actually returned.""" + stats = results.get("stats") or {} + solve = results.get("solve") + reference = results.get("reference") + palette = session.palette + now = datetime.now(timezone.utc).strftime("%Y-%m-%d %H:%M UTC") + + lines = [] + add = lines.append + + add(f"# {session.target} - {session.telescope}") + add("") + add(f"Processed {now} by the pipeline in the " + "[astrophotography](https://git.discworld.casa/laurence/astrophotography) " + "repository (`pipeline/`). This file is generated: re-running the " + "pipeline rewrites it from the actual results.") + add("") + + # ---------------------------------------------------------- what is here + add("## What was captured") + add("") + add("| Filter | Frames | Exposure | Binning | Total |") + add("|---|---|---|---|---|") + total = 0.0 + for filt in session.filters: + frames = session.by_filter(filt) + exps = sorted({f.exptime for f in frames}) + bins = sorted({f.binning for f in frames}) + secs = sum(f.exptime for f in frames) + total += secs + add(f"| {filt} | {len(frames)} | " + f"{'/'.join(f'{e:.0f} s' for e in exps)} | " + f"{'/'.join(f'bin{b}' for b in bins)} | {_fmt_exposure(secs)} |") + add(f"| **Total** | **{len(session.frames)}** | | | " + f"**{_fmt_exposure(total)}** |") + add("") + add(f"Colour scheme: **{palette}**" + + {"SHO": " - the Hubble palette, SII to red, Ha to green, OIII to blue", + "HOO": " - Ha to red, OIII to green and blue", + "LRGB": " - colour from R/G/B, detail from the deeper luminance", + "RGB": " - no luminance channel, so brightness is synthesised from " + "the colour channels", + "MONO": " - a single filter, rendered as greyscale"}.get(palette, "")) + add("") + + # -------------------------------------------------------- frame quality + if stats: + add("## Frame quality") + add("") + add("| Filter | Sky (ADU) | Noise (ADU) | Seeing (px) |") + add("|---|---|---|---|") + for filt in session.filters: + vals = [stats[f.name] for f in session.by_filter(filt) + if f.name in stats] + if not vals: + continue + sky = np.median([v["sky"] for v in vals]) + rms = np.median([v["rms"] for v in vals]) + fwhm = np.median([v["fwhm"] for v in vals]) + add(f"| {filt} | {sky:.1f} | {rms:.2f} | {fwhm:.2f} |") + add("") + + # ---------------------------------------------------------- what was done + add("## What was done") + add("") + add("1. **Ingest.** Archives uncompressed into `calibrated/`. Only the " + "`calibrated-` set is extracted; the `raw-` archives stay zipped in " + "`raw/` as untouched originals, because extracting both would put two " + "copies of every frame into the stack.") + add("2. **Measure.** Sky, noise and star sharpness per frame, with each " + "frame's star list cached for the stages that follow.") + if reference: + add(f"3. **Register and stack.** Every frame aligned by matching star " + f"patterns onto one reference - `{reference.name}` " + f"({reference.filter}, the sharpest frame of the deepest filter) - " + f"so all masters share a pixel grid and the colour image needs no " + f"further alignment. Combined with a " + f"{'sigma-clipped ' if len(session.frames) >= 3 else ''}" + f"noise-weighted mean.") + else: + add("3. **Register and stack.**") + if solve: + add(f"4. **Plate solve.** Matched against Gaia DR3: " + f"**{solve['nstars']} stars, {solve['residual_arcsec']:.2f} arcsec " + f"residual**. The WCS is written into every master, so any position " + f"measured from these images is quotable.") + else: + add("4. **Plate solve. FAILED for this session** - the star pattern " + "match did not converge against Gaia. The images are unaffected, " + "but nothing here carries sky coordinates, so positions cannot be " + "measured from them.") + add(f"5. **Colour.** Assembled as {palette}, with the sky background " + "levelled by a plane fitted outside the target and the three channels' " + "backgrounds forced neutral.") + add("") + + # ------------------------------------------------------------ astrometry + if solve: + add("## Astrometry") + add("") + add(f"- Field centre: **{solve['centre']}**") + add(f"- Plate scale: **{solve['scale']:.4f} arcsec/pixel**") + add(f"- Field of view: {solve['fov_arcmin'][0]:.1f}' x " + f"{solve['fov_arcmin'][1]:.1f}'") + add(f"- Position angle: {solve['position_angle']:.2f} deg") + add(f"- Residual: {solve['residual_px']:.2f} px " + f"({solve['residual_arcsec']:.2f} arcsec) on {solve['nstars']} stars") + add("") + + # ----------------------------------------------------------- honest limits + add("## What limits this data") + add("") + shortest = min((sum(f.exptime for f in session.by_filter(x)) + for x in session.filters), default=0) + if len(session.frames) < 6: + add(f"- **Very few frames ({len(session.frames)} in total).** With so " + "few there is no meaningful outlier rejection, so satellite " + "trails and cosmic rays survive into the masters.") + if shortest < 600: + add(f"- **Short integration** - the thinnest channel has only " + f"{_fmt_exposure(shortest)}. No processing recovers signal that " + "was never collected.") + if palette == "RGB": + add("- **No luminance channel**, so brightness is synthesised from the " + "colour data and the result is noisier than an LRGB set of the " + "same total exposure would be.") + if palette == "MONO": + add("- **A single filter**, so there is no colour information at all.") + if not solve: + add("- **No astrometric solution**, so nothing here can be used for " + "positions or cross-matched against a catalogue.") + fwhms = [v["fwhm"] for v in stats.values() if np.isfinite(v.get("fwhm", np.nan))] + if fwhms and solve: + seeing = float(np.median(fwhms)) * solve["scale"] + add(f"- **Seeing of about {seeing:.1f} arcsec**, which sets the finest " + "detail present regardless of how the data is processed.") + add("") + + # ------------------------------------------------------------- directories + add("## Where things are") + add("") + add("| Directory | Contents |") + add("|---|---|") + add("| `raw/` | exactly what the telescope delivered; untouched |") + add("| `calibrated/` | the uncompressed calibrated frames |") + add("| `stacks/masters/` | one registered, combined master per filter |") + add("| `final/` | the finished image |") + add("| `intermediates/` | caches; deletable, regenerated by a re-run |") + add("") + add("Re-run any stage by pointing the pipeline at this directory:") + add("") + add("```") + add(f"set ASTRO_SESSION={session.root}") + add("python pipeline/run.py") + add("```") + add("") + + path = os.path.join(session.root, "METHODS.md") + with open(path, "w", encoding="utf-8") as fh: + fh.write("\n".join(lines)) + if verbose: + print(f" -> METHODS.md") + return path diff --git a/pipeline/run.py b/pipeline/run.py new file mode 100644 index 0000000..e58cc07 --- /dev/null +++ b/pipeline/run.py @@ -0,0 +1,158 @@ +"""The push-button entry point: process a session, or a whole archive. + + python run.py # process $ASTRO_SESSION + python run.py PATH # process one session + python run.py --all PATH # find and process every session + python run.py --stage register PATH # run one stage + python run.py --force PATH # ignore existing outputs + +Stages run in order and each is skipped when its output already exists, so an +interrupted run continues rather than starting again. A failure in one stage of +one session does not stop the others: the failure is recorded and the batch +carries on, because in a batch of twenty the useful outcome is nineteen results +and one clear error, not nothing. +""" +import argparse +import os +import sys +import time +import traceback + +import layout +import session as session_mod + +STAGES = ("ingest", "measure", "register", "solve", "colour", "report") + + +def _exists(path): + return os.path.exists(path) and os.path.getsize(path) > 0 + + +def process(root, stages=STAGES, force=False, verbose=True): + """Run the requested stages for one session. Returns a result summary.""" + s = session_mod.Session(root) + result = {"root": root, "stages": {}, "ok": True} + t0 = time.time() + + if "ingest" in stages: + s.ingest(verbose=verbose) + s.scan(verbose=verbose) + if not s.frames: + print(" no light frames found - skipping") + result["ok"] = False + result["error"] = "no frames" + return result + s.save_manifest() + result["target"] = s.target + result["telescope"] = s.telescope + result["palette"] = s.palette + result["frames"] = len(s.frames) + + results = {} + try: + if "measure" in stages: + import measure + cache = os.path.join(root, layout.INTERMEDIATES, "_measure.npz") + if force or not _exists(cache): + print(" measure") + measure.run(s, force=force, verbose=False) + else: + print(" measure (cached)") + _, results["stats"] = measure.load(s) + + if "register" in stages: + import measure + import register + masters = os.path.join(root, layout.MASTERS) + have = _exists(os.path.join(masters, f"master-{s.filters[0]}.fit")) + if force or not have: + print(" register") + _, ref = register.run(s, verbose=verbose) + else: + print(" register (cached)") + _, st = measure.load(s) + ref = measure.choose_reference(s, st) + results["reference"] = ref + + if "solve" in stages: + import astrometry + print(" solve") + try: + results["solve"] = astrometry.run(s, verbose=verbose) + except Exception as exc: # noqa: BLE001 + # An unsolved session still produces images; record and go on. + print(f" solve failed: {type(exc).__name__}: {exc}") + results["solve"] = None + + if "colour" in stages: + import colour + print(" colour") + results["final"] = colour.run(s, verbose=verbose) + + if "report" in stages: + import report + print(" report") + report.write(s, results, verbose=verbose) + + except Exception as exc: # noqa: BLE001 + result["ok"] = False + result["error"] = f"{type(exc).__name__}: {exc}" + print(f" FAILED: {result['error']}") + if verbose: + traceback.print_exc(limit=3) + + result["solved"] = bool(results.get("solve")) + result["seconds"] = time.time() - t0 + return result + + +def main(argv=None): + ap = argparse.ArgumentParser(description=__doc__) + ap.add_argument("path", nargs="?", default=None, + help="session directory, or a tree with --all") + ap.add_argument("--all", action="store_true", + help="find and process every session beneath PATH") + ap.add_argument("--stage", action="append", choices=STAGES, + help="run only this stage (repeatable)") + ap.add_argument("--force", action="store_true", + help="re-run stages even if their output exists") + args = ap.parse_args(argv) + + root = args.path or layout.SESSION + stages = tuple(args.stage) if args.stage else STAGES + + targets = session_mod.discover(root) if args.all else [root] + if not targets: + print(f"no sessions found under {root}") + return 1 + print(f"{len(targets)} session(s), stages: {', '.join(stages)}\n") + + results = [] + for i, t in enumerate(targets, 1): + print(f"[{i}/{len(targets)}] {os.path.relpath(t, root)}") + results.append(process(t, stages=stages, force=args.force)) + print() + + print("=" * 72) + print(f"{'target':16s} {'scope':16s} {'palette':8s} " + f"{'frames':>6s} {'solved':>7s} {'time':>7s}") + for r in results: + if not r.get("ok") and "target" not in r: + print(f"{os.path.basename(r['root']):16s} " + f"FAILED: {r.get('error', 'unknown')}") + continue + print(f"{r.get('target', '?')[:16]:16s} " + f"{r.get('telescope', '?')[:16]:16s} " + f"{r.get('palette', '?'):8s} {r.get('frames', 0):6d} " + f"{'yes' if r.get('solved') else 'no':>7s} " + f"{r.get('seconds', 0):6.0f}s") + failed = [r for r in results if not r.get("ok")] + unsolved = [r for r in results if r.get("ok") and not r.get("solved")] + print(f"\n{len(results) - len(failed)}/{len(results)} processed" + + (f", {len(unsolved)} without an astrometric solution" if unsolved + else "")) + return 1 if failed else 0 + + +if __name__ == "__main__": + sys.exit(main())