panoptes.pipeline.utils.images¶

def extract_metadata(header, path_info, camera_settings: CameraSettings | None = None) -> dict:
    """Get the metadata from a FITS image.

    Returns ``dict(unit=..., sequence=..., image=...)``, one key per document.
    That shape is the contract and predates this function's current form; see
    data contract 3.1.

    ``camera_settings`` supplies the fallbacks for calibration values the
    header does not carry. Whatever it supplies is labeled `Provenance.DEFAULT`
    in the resulting document rather than being indistinguishable from a
    measurement.
    """
    camera_settings = camera_settings or CameraSettings()
    try:
        measured_rggb = [float(x) for x in header.get("MEASRGGB", "0 0 0 0").split(" ")]
        calibration = provenance.resolve_camera(header, camera_settings)

        # `str(header.get("CAMSN"))` produced the literal string "None" on every
        # frame POCS wrote without the keyword, which reads downstream as a body
        # serial that 37 cameras share. Absent is absent.
        camera_serial = header.get("CAMSN")
        camera_serial = str(camera_serial) if camera_serial is not None else None
        file_date = path_info.image_time.to_datetime(timezone=UTC)
        camera_date = parse_date(header.get("DATE-OBS", path_info.image_time)).replace(tzinfo=UTC)

        unit_info = dict(
            unit_id=path_info.unit_id,
            latitude=header.get("LAT-OBS"),
            longitude=header.get("LONG-OBS"),
            elevation=float(header.get("ELEV-OBS")),
            name=header.get("OBSERVER"),
        )

        sequence_info = dict(
            sequence_id=path_info.sequence_id,
            sequence_time=path_info.sequence_time.to_datetime(timezone=UTC),
            coordinates=dict(
                airmass=header.get("AIRMASS"),
                mount_dec=header.get("DEC-MNT"),
                mount_ra=header.get("RA-MNT"),
                mount_ha=header.get("HA-MNT"),
            ),
            camera=dict(
                camera_id=path_info.camera_id,
                lens_serial_number=header.get("INTSN"),
                serial_number=camera_serial,
            ),
            # NAXIS1/NAXIS2, not IMAGEW/IMAGEH: see `provenance.resolve_camera`.
            imagew=calibration["image_width"].value,
            imageh=calibration["image_height"].value,
            field_name=header.get("FIELD", ""),
            software_version=header.get("CREATOR", ""),
        )

        image_info = dict(
            uid=path_info.get_full_id(sep="_"),
            camera=dict(
                # The camera identifiers repeat here on purpose. They sat only
                # in `sequence_info`, which made joining a per-frame white level
                # to a per-camera profile awkward -- the measurement is per
                # frame and the identifier was per sequence. Both are in every
                # header, so carrying them costs nothing. See data contract 3.1.
                camera_id=path_info.camera_id,
                serial_number=camera_serial,
                blue_balance=float(header.get("BLUEBAL")),
                circconf=float(header.get("CIRCCONF", "0.").split(" ")[0]),
                colortemp=float(header.get("COLORTMP")),
                dateobs=camera_date,
                exptime=float(header.get("EXPTIME")),
                iso=header.get("ISO"),
                measured_b=measured_rggb[3],
                measured_ev1=float(header.get("MEASEV")),
                measured_ev2=float(header.get("MEASEV2")),
                measured_g1=measured_rggb[1],
                measured_g2=measured_rggb[2],
                measured_r=measured_rggb[0],
                red_balance=float(header.get("REDBAL")),
                temperature=float(header.get("CAMTEMP", 0).split(" ")[0]),
                white_lvln=header.get("WHTLVLN"),
                white_lvls=header.get("WHTLVLS"),
            ),
            environment=dict(
                moonfrac=float(header.get("MOONFRAC")),
                moonsep=float(header.get("MOONSEP")),
            ),
            # Resolved calibration, each value with the tier it came from, so a
            # header reading and a fleet-wide default do not serialize alike.
            calibration=provenance.as_document(calibration),
            file_creation_date=file_date,
            image_time=path_info.image_time.to_datetime(timezone=UTC),
        )

    except Exception as e:
        print(f"Error in extracting metadata: {e!r}")
        raise e

    if fell_back := provenance.defaulted(calibration):
        logger.warning(
            f"{path_info.id}: no better source than a fleet-wide default for "
            f"{', '.join(fell_back)} -- see improvement plan 1.4"
        )

    print("Metadata extracted from header")
    return dict(unit=unit_info, sequence=sequence_info, image=image_info)

def match_sources(
    detected_sources: pandas.DataFrame,
    solved_wcs0: WCS,
    settings: ImageSettings,
    image_edge: int = 10,
    image_width: int | None = None,
    image_height: int | None = None,
) -> pandas.DataFrame:
    """Match detected sources against the catalog and drop the ones near the edge.

    `image_width` and `image_height` are passed in rather than read from
    `settings` because the settings hold fleet-wide defaults, and the frame's
    real dimensions come from its own `NAXIS1`/`NAXIS2`. The old notebook
    reconciled this by writing the shape back into the settings object, which
    would now change the params fingerprint from frame to frame and make the
    idempotency check meaningless.
    """
    image_width = image_width or settings.params.camera.image_width
    image_height = image_height or settings.params.camera.image_height
    print(f"Matching {len(detected_sources)} sources to wcs.")
    catalog_filename = settings.params.catalog.catalog_filename
    vmag_limits = settings.params.catalog.vmag_limits
    print(f"Using catalog from {catalog_filename}")
    catalog_sources = sources.get_stars_from_wcs(
        solved_wcs0,
        catalog_filename=catalog_filename,
        vmag_min=vmag_limits[0],
        vmag_max=vmag_limits[1],
    )
    print(f"Matching sources to catalog for {len(detected_sources)} sources")
    matched_sources = sources.get_catalog_match(
        detected_sources,
        wcs=solved_wcs0,
        catalog_stars=catalog_sources,
        ra_column="photutils_sky_centroid_ra",
        dec_column="photutils_sky_centroid_dec",
        max_separation_arcsec=settings.params.catalog.max_separation_arcsec,
    )
    # Drop matches near border
    print(f"Filtering sources near within {image_edge} pixels of {image_width}x{image_height}")
    matched_sources = matched_sources.query(
        f"catalog_wcs_x_int > {image_edge} and "
        f"catalog_wcs_x_int < {image_width - image_edge} and "
        f"catalog_wcs_y_int > {image_edge} and "
        f"catalog_wcs_y_int < {image_height - image_edge}"
    ).copy()
    print(f"Found {len(matched_sources)} matching sources")

    # There should not be too many duplicates at this point and they are returned in order
    # of catalog separation, so we take the first.
    duplicates = matched_sources.duplicated("picid", keep="first")
    print(f"Found {len(matched_sources[duplicates])} duplicate sources")

    # Mark which ones were duplicated.
    matched_sources.loc[:, "catalog_match_duplicate"] = False
    dupes = matched_sources.picid.isin(matched_sources[duplicates].picid)
    matched_sources.loc[dupes, "catalog_match_duplicate"] = True

    # Filter out duplicates.
    matched_sources = matched_sources.loc[~duplicates].copy()
    print(f"Found {len(matched_sources)} matching sources after removing duplicates")

    # Add some binned fields that we can use for table partitioning.
    matched_sources["catalog_dec_bin"] = matched_sources.catalog_dec.astype("int")
    matched_sources["catalog_ra_bin"] = matched_sources.catalog_ra.astype("int")

    # Color excesses, where the catalog carries Gaia photometry. These are an
    # enrichment, not a requirement: a catalog holding only the four columns
    # `REQUIRED_CATALOG_COLUMNS` names is valid, and AGENTS.md promises it can
    # be rebuilt without touching this package. Computing them unconditionally
    # broke that promise -- such a catalog passed `read_catalog`, matched, and
    # then died on `catalog_gaiabp` after solving and detection had already run.
    #
    # Nothing reads these columns yet. They are here for reference selection,
    # where matching a reference's color to the target's is the point, and that
    # step does not exist. See #202.
    if set(sources.GAIA_CATALOG_COLUMNS) <= set(matched_sources.columns):
        matched_sources["catalog_gaia_bg_excess"] = (
            matched_sources.catalog_gaiabp - matched_sources.catalog_gaiamag
        )
        matched_sources["catalog_gaia_br_excess"] = (
            matched_sources.catalog_gaiabp - matched_sources.catalog_gaiarp
        )
        matched_sources["catalog_gaia_rg_excess"] = (
            matched_sources.catalog_gaiarp - matched_sources.catalog_gaiamag
        )
    else:
        logger.debug(
            "Catalog carries no Gaia photometry; color excesses are not computed. "
            f"Add {list(sources.GAIA_CATALOG_COLUMNS)} to the catalog to get them."
        )

    return matched_sources

def plate_solve(settings: ImageSettings, filename=None, timeout=30, replace=False, **kwargs):
    """Solve `filename` and return its WCS.

    `replace` defaults to False, which is `panoptes-utils`' non-destructive
    mode: the WCS goes to a `.new` file and the input is left alone. The
    library default is True, and it does not merely add a WCS -- `solve-field`
    rewrites the file as a single HDU, so solving a multi-extension product in
    place silently discards every extension beside the science frame. Verified
    against a real frame: a four-extension file came back with one.

    The same applies to `panoptes-utils image fits solve`, which forwards its
    keywords to the same function and so carries the same default.
    """
    filename = filename or settings.files.reduced_filename
    print(f"Plate solving {filename}")

    # Add custom options for solving.
    options = [
        "--fits-image",
        "--scale-low",
        "10",
        "--scale-high",
        "20",
        "--scale-units",
        "degw",
        "--radius",
        "15",
        "--guess-scale",
        "--no-background-subtraction",
        "--cpulimit",
        str(timeout),
        "--no-verify",
        "--crpix-center",
        "--temp-axy",
        "--index-xyls",
        "none",
        "--solved",
        "none",
        "--match",
        "none",
        "--rdls",
        "none",
        "--corr",
        filename.with_suffix(".corr").as_posix(),
        "--downsample",
        "4",
        "--tweak-order",
        "4",
        "--no-plots",
    ]

    solved_headers = fits_utils.get_solve_field(
        str(filename),
        skip_solved=False,
        replace=replace,
        solve_opts=options,
        timeout=timeout,
        **kwargs,
    )
    solved_path = solved_headers.pop("solved_fits_file")
    print(f"Solving completed successfully for {solved_path}")
    solved_wcs0 = WCS(solved_headers)
    return solved_wcs0