search

neuromorphicsystems / astrometry

Astrometry turns a list of star positions into a pixel-to-sky transformation (WCS)
GNU General Public License v3.0
27 stars 3 forks source link

readme

Astrometry

Astrometry turns a list of star positions into a pixel-to-sky transformation (WCS) by calling C functions from the Astrometry.net library (https://astrometry.net).

Astrometry.net star index files ("series") are automatically downloaded when required.

This package is useful for solving plates from a Python script, comparing star extraction methods, or hosting a simple local version of Astrometry.net with minimal dependencies. See https://github.com/dam90/astrometry for a more complete self-hosting solution.

Unlike Astrometry.net, Astrometry does not include FITS parsing or image pre-processing algorithms. Stars must be provided as a list of pixel positions.

This library works on Linux and macOS, but not Windows (at the moment). WSL should work but has not been tested.

We are not the authors of the Astrometry.net library. You should cite works from https://astrometry.net/biblio.html if you use the Astrometry.net algorithm via this package.

Get started

python3 -m pip install astrometry
import astrometry

solver = astrometry.Solver(
    astrometry.series_5200.index_files(
        cache_directory="astrometry_cache",
        scales={6},
    )
)

stars = [
    [388.9140568247906, 656.5003281719216],
    [732.9210858972549, 473.66395545775106],
    [401.03459504299843, 253.788113189415],
    [312.6591868096163, 624.7527729425295],
    [694.6844564647456, 606.8371776658344],
    [741.7233477959561, 344.41284826261443],
    [867.3574610200455, 672.014835980283],
    [1063.546651153479, 593.7844603550848],
    [286.69070190952704, 422.170016812049],
    [401.12779619355155, 16.13543616977013],
    [205.12103484692776, 698.1847350789413],
    [202.88444768690894, 111.24830187635557],
    [339.1627757703069, 86.60739435924549],
]

solution = solver.solve(
    stars=stars,
    size_hint=None,
    position_hint=None,
    solution_parameters=astrometry.SolutionParameters(),
)

if solution.has_match():
    print(f"{solution.best_match().center_ra_deg=}")
    print(f"{solution.best_match().center_dec_deg=}")
    print(f"{solution.best_match().scale_arcsec_per_pixel=}")

solve is thread-safe. It can be called any number of times from the same Solver object.

Examples

Provide size and position hints

import astrometry

solver = ...
solution = solver.solve(
    stars=...,
    size_hint=astrometry.SizeHint(
        lower_arcsec_per_pixel=1.0,
        upper_arcsec_per_pixel=2.0,
    ),
    position_hint=astrometry.PositionHint(
        ra_deg=65.7,
        dec_deg=36.2,
        radius_deg=1.0,
    ),
    solution_parameters=...
)

Print progress information (download and solve)

import astrometry
import logging

logging.getLogger().setLevel(logging.INFO)

solver = ...
solution = ...

Print field stars metadata

Astrometry extracts metadata from the star index ("series"). See Choosing series for a description of the available data.

import astrometry

solver = ...
solution = ...

if solution.has_match():
    for star in solution.best_match().stars:
        print(f"{star.ra_deg}º, {star.dec_deg}º:", star.metadata)

Calculate field stars pixel positions with astropy

import astrometry

solver = ...
solution = ...

if solution.has_match():
    wcs = solution.best_match().astropy_wcs()
    pixels = wcs.all_world2pix(
        [[star.ra_deg, star.dec_deg] for star in solution.best_match().stars],
        0,
    )
    # pixels is a len(solution.best_match().stars) x 2 numpy array of float values

astropy.wcs.WCS provides many more functions to probe the transformation properties and convert from and to pixel coordinates. See https://docs.astropy.org/en/stable/api/astropy.wcs.WCS.html for details. Astropy (https://pypi.org/project/astropy/) must be installed to use this method.

Print series description and size (without downloading them)

import astrometry

print(astrometry.series_5200_heavy.description)
print(astrometry.series_5200_heavy.size_as_string({2, 3, 4}))

See Choosing Series for a list of available series.

Disable tune-up and distortion

import astrometry

solver = ...
solution = solver.solve(
    stars=...,
    size_hint=...,
    position_hint=...,
    solution_parameters=astrometry.SolutionParameters(
        sip_order=0,
        tune_up_logodds_threshold=None,
    ),
)

Stop the solver early using the log-odds callback

Return after the first match

import astrometry

solver = ...
solution = solver.solve(
    stars=...,
    size_hint=...,
    position_hint=...,
    solution_parameters=astrometry.SolutionParameters(
        logodds_callback=lambda logodds_list: astrometry.Action.STOP,
    ),
)

Return early if the best log-odds are larger than 100.0

import astrometry

solver = ...
solution = solver.solve(
    stars=...,
    size_hint=...,
    position_hint=...,
    solution_parameters=astrometry.SolutionParameters(
        logodds_callback=lambda logodds_list: (
            astrometry.Action.STOP
            if logodds_list[0] > 100.0
            else astrometry.Action.CONTINUE
        ),
    ),
)

Return early if there are at least ten matches

import astrometry

solver = ...
solution = solver.solve(
    stars=...,
    size_hint=...,
    position_hint=...,
    solution_parameters=astrometry.SolutionParameters(
        logodds_callback=lambda logodds_list: (
            astrometry.Action.STOP
            if len(logodds_list) >= 10.0
            else astrometry.Action.CONTINUE
        ),
    ),
)

Return early if the three best matches are similar

import astrometry

def logodds_callback(logodds_list: list[float]) -> astrometry.Action:
    if len(logodds_list) < 3:
        return astrometry.Action.CONTINUE
    if logodds_list[1] > logodds_list[0] - 10 and logodds_list[2] > logodds_list[0] - 10:
        return astrometry.Action.STOP
    return astrometry.Action.CONTINUE

solver = ...
solution = solver.solve(
    stars=...,
    size_hint=...,
    position_hint=...,
    solution_parameters=astrometry.SolutionParameters(
        logodds_callback=logodds_callback,
    ),
)

Choosing series

This library downloads series from http://data.astrometry.net. A solver can be instantiated with multiple series and scales as follows:

import astrometry

solver = astrometry.Solver(
    astrometry.series_5200.index_files(
        cache_directory="astrometry_cache",
        scales={4, 5, 6},
    )
    + astrometry.series_4200.index_files(
        cache_directory="astrometry_cache",
        scales={6, 7, 12},
    )
)

Astrometry.net gives the following recommendations to choose a scale:

Each index file contains a large number of “skymarks” (landmarks for the sky) that allow our solver to identify your images. The skymarks contained in each index file have sizes (diameters) within a narrow range. You probably want to download index files whose quads are, say, 10% to 100% of the sizes of the images you want to solve.

For example, let’s say you have some 1-degree square images. You should grab index files that contain skymarks of size 0.1 to 1 degree, or 6 to 60 arcminutes. Referring to the table below, you should [try index files with scales 3 to 9]. You might find that the same number of fields solve, and faster, using just one or two of the index files in the middle of that range - in our example you might try [5, 6 and 7].

-- http://astrometry.net/doc/readme.html

Scale 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19
Skymark diameter (arcmin) [2.0, 2.8] [2.8, 4.0] [4.0, 5.6] [5.6, 8.0] [8, 11] [11, 16] [16, 22] [22, 30] [30, 42] [42, 60] [60, 85] [85, 120] [120, 170] [170, 240] [240, 340] [340, 480] [480, 680] [680, 1000] [1000, 1400] [1400, 2000]

The table below lists series sizes and properties (we copied the descriptions from http://data.astrometry.net). You can access a series' object with astrometry.series_{name} (for example astrometry.series_4200).

Name Total size Scales Description Metadata
4100 0.36 GB [7, 19] built from the Tycho-2 catalog, good for images wider than 1 degree, recommended MAG_BT: float
MAG_VT: float
MAG_HP: float
MAG: float
4200 33.78 GB [0, 19] built from the near-infared 2MASS survey, runs out of stars at the low end, most users will probably prefer 4100 or 5200 j_mag: float
5000 76.24 GB [0, 7] an older version from Gaia-DR2 but without Tycho-2 stars merged in, our belief is that series_5200 will work better than this one source_id: int
phot_g_mean_mag: float
phot_bp_mean_mag: float
phot_rp_mean_mag: float
parallax: float
parallax_error: float
pmra: float
pmra_error: float
pmdec: float
pmdec_error: float
ra: float
dec: float
ref_epoch: float
5200 36.14 GB [0, 6] LIGHT version built from Tycho-2 + Gaia-DR2, good for images narrower than 1 degree, combine with 4100-series for broader scale coverage, the LIGHT version contains smaller files with no additional Gaia-DR2 information tagged along, recommended -
5200_heavy 79.67 GB [0, 6] HEAVY version same as 5200, but with additional Gaia-DR2 information (magnitude in G, BP, RP, proper motions and parallaxes), handy if you want that extra Gaia information for matched stars ra: float
dec: float
mag: float
ref_cat: str
ref_id: int
pmra: float
pmdec: float
parallax: float
ra_ivar: float
dec_ivar: float
pmra_ivar: float
pmdec_ivar: float
parallax_ivar: float
phot_bp_mean_mag: float
phot_rp_mean_mag: float
6000 1.20 GB [4, 6] very specialized, uses GALEX Near-UV measurements, and only a narrow range of scales fuv_mag: float
nuv_mag: float
6100 1.58 GB [4, 6] very specialized, uses GALEX Far-UV measurements, and only a narrow range of scales fuv_mag: float
nuv_mag: float

The table below indicates the total file size for each scale (most series have multiple index files per scale).

Name 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19
4100 - - - - - - - 165.00 MB 94.55 MB 49.77 MB 24.87 MB 10.21 MB 5.30 MB 2.73 MB 1.38 MB 740.16 kB 408.96 kB 247.68 kB 187.20 kB 144.00 kB
4200 14.22 GB 9.25 GB 5.06 GB 2.63 GB 1.31 GB 659.09 MB 328.25 MB 165.44 MB 81.84 MB 41.18 MB 20.52 MB 8.02 MB 4.17 MB 2.16 MB 1.10 MB 596.16 kB 339.84 kB 213.12 kB 164.16 kB 132.48 kB
5000 34.79 GB 20.19 GB 10.74 GB 5.44 GB 2.71 GB 1.36 GB 676.79 MB 340.73 MB - - - - - - - - - - - -
5200 17.20 GB 9.49 GB 4.86 GB 2.45 GB 1.22 GB 614.89 MB 307.72 MB - - - - - - - - - - - - -
5200_heavy 36.46 GB 21.20 GB 11.29 GB 5.72 GB 2.85 GB 1.43 GB 714.56 MB - - - - - - - - - - - - -
6000 - - - - 892.55 MB 457.66 MB 233.23 MB - - - - - - - - - - - - -
6100 - - - - 599.33 MB 384.09 MB 214.79 MB - - - - - - - - - - - - -

Documentation

Solver

class Solver:
    def __init__(self, index_files: list[pathlib.Path]): ...

    def solve(
        self,
        stars: typing.Iterable[SupportsFloatMapping],
        size_hint: typing.Optional[SizeHint],
        position_hint: typing.Optional[PositionHint],
        solution_parameters: SolutionParameters,
    ) -> Solution: ...

solve is thread-safe and can be called any number of times.

SizeHint

@dataclasses.dataclass
class SizeHint:
    lower_arcsec_per_pixel: float
    upper_arcsec_per_pixel: float

lower_arcsec_per_pixel and upper_arcsec_per_pixel must be larger than 0 and upper_arcsec_per_pixel must be smaller than or equal to upper_arcsec_per_pixel.

PositionHint

@dataclasses.dataclass
class PositionHint:
    ra_deg: float
    dec_deg: float
    radius_deg: float

All values are in degrees and must use the same frame of reference as the index files. Astrometry.net index files use J2000 FK5 (https://docs.astropy.org/en/stable/api/astropy.coordinates.FK5.html). ICRS and FK5 differ by less than 0.1 arcsec (https://www.iers.org/IERS/EN/Science/ICRS/ICRS.html).

Action

class Action(enum.Enum):
    STOP = 0
    CONTINUE = 1

Parity

class Parity(enum.IntEnum):
    NORMAL = 0
    FLIP = 1
    BOTH = 2

SolutionParameters

@dataclasses.dataclass
class SolutionParameters:
    solve_id: typing.Optional[str] = None
    uniformize_index: bool = True
    deduplicate: bool = True
    sip_order: int = 3
    sip_inverse_order: int = 0
    distance_from_quad_bonus: bool = True
    positional_noise_pixels: float = 1.0
    distractor_ratio: float = 0.25
    code_tolerance_l2_distance: float = 0.01
    minimum_quad_size_pixels: typing.Optional[float] = None
    minimum_quad_size_fraction: float = 0.1
    maximum_quad_size_pixels: float = 0.0
    maximum_quads: int = 0
    maximum_matches: int = 0
    parity: Parity = Parity.BOTH
    tune_up_logodds_threshold: typing.Optional[float] = 14.0
    output_logodds_threshold: float = 21.0
    slices_generator: typing.Callable[[int], typing.Iterable[tuple[int, int]]] = astrometry.batches_generator(25)
    logodds_callback: typing.Callable[[list[float]], Action] = lambda _: Action.CONTINUE

Accepted matches are always tuned up, even if they hit tune_up_logodds_threshold and were already tuned-up. Since log-odds are compared with the thresholds before the tune-up, the final log-odds are often significantly larger than output_logodds_threshold. Set tune_up_logodds_threshold to a value larger than or equal to output_logodds_threshold to disable the first tune-up, and None to disable tune-up altogether. Tune-up logic is equivalent to the following Python snippet:

# This (pseudo-code) snippet assumes the following definitions:
# match: candidate match object
# log_odds: current match log-odds
# add_to_solution: appends the match to the solution list
# tune_up: tunes up a match object and returns the new match and the new log-odds
if tune_up_logodds_threshold is None:
    if log_odds >= output_logodds_threshold:
        add_to_solution(match)
else:
    if log_odds >= output_logodds_threshold:
        tuned_up_match, tuned_up_loggods = tune_up(match)
        add_to_solution(tuned_up_match)
    elif log_odds >= tune_up_logodds_threshold:
        tuned_up_match, tuned_up_loggods = tune_up(match)
        if tuned_up_loggods >= output_logodds_threshold:
            tuned_up_twice_match, tuned_up_twice_loggods = tune_up(tuned_up_match)
            add_to_solution(tuned_up_twice_match)

Astrometry.net gives the following description of the tune-up algorithm. See tweak2 in astrometry.net/solver/tweak2.c for the implementation.

Given an initial WCS solution, compute SIP polynomial distortions using an annealing-like strategy. That is, it finds matches between image and reference catalog by searching within a radius, and that radius is small near a region of confidence, and grows as you move away. That makes it possible to pick up more distant matches, but they are downweighted in the fit. The annealing process reduces the slope of the growth of the matching radius with respect to the distance from the region of confidence.

-- astrometry.net/include/astrometry/tweak2.h

Solution

@dataclasses.dataclass
class Solution:
    solve_id: str
    matches: list[Match]

    def has_match(self) -> bool: ...

    def best_match(self) -> Match: ...

    def to_json(self) -> str: ...

    @classmethod
    def from_json(cls, solution_as_json: str) -> Solution: ...

matches are sorted in descending log-odds order. best_match returns the first match in the list. to_json and from_json may be used to save and load solutions.

Match

@dataclasses.dataclass
class Match:
    logodds: float
    center_ra_deg: float
    center_dec_deg: float
    scale_arcsec_per_pixel: float
    index_path: pathlib.Path
    stars: tuple[Star, ...]
    quad_stars: tuple[Star, ...]
    wcs_fields: dict[str, tuple[typing.Any, str]]

    def astropy_wcs(self) -> astropy.wcs.WCS: ...

astropy_wcs generates an Astropy WCS object. Astropy (https://pypi.org/project/astropy/) must be installed to use this method. See Calculate field stars pixel positions with astropy for details.

Star

@dataclasses.dataclass
class Star:
    ra_deg: float
    dec_deg: float
    metadata: dict[str, typing.Any]

ra_deg and dec_deg are in degrees and use the same frame of reference as the index files. Astrometry.net index files use J2000 FK5 (https://docs.astropy.org/en/stable/api/astropy.coordinates.FK5.html). ICRS and FK5 differ by less than 0.1 arcsec (https://www.iers.org/IERS/EN/Science/ICRS/ICRS.html).

The contents of metadata depend on the data available in index files. See Series for details.

Series

@dataclasses.dataclass
class Series:
    name: str
    description: str
    scale_to_sizes: dict[int, tuple[int, ...]]
    url_pattern: str

    def size(self, scales: typing.Optional[set[int]] = None): ...

    def size_as_string(self, scales: typing.Optional[set[int]] = None): ...

    def index_files(
        self,
        cache_directory: typing.Union[bytes, str, os.PathLike],
        scales: typing.Optional[set[int]] = None,
    ) -> list[pathlib.Path]: ...

Change the constants astrometry.CHUNK_SIZE, astrometry.DOWNLOAD_SUFFIX and astrometry.TIMEOUT to configure the downloader parameters.

batches_generator

def batches_generator(
    batch_size: int,
) -> typing.Callable[[int], typing.Iterable[tuple[int, int]]]:
    ...

batches_generator returns a slices generator compatible with SolutionParameters.slices_generator. The slices are non-overlapping and non-full slices are ignored. For instance, a batch size of 25 over 83 stars would generate the slices (0, 25), (25, 50), and (50, 75).

SupportsFloatMapping

class SupportsFloatMapping(typing.Protocol):
    def __getitem__(self, index: typing.SupportsIndex, /) -> typing.SupportsFloat:
        ...

Contribute

Clone this repository and pull its submodule:

git clone --recursive https://github.com/neuromorphicsystems/astrometry.git
cd astrometry

or

git clone https://github.com/neuromorphicsystems/astrometry.git
cd astrometry
git submodule update --recursive

Format the code:

clang-format -i astrometry_extension/astrometry_extension.c astrometry_extension/astrometry_extension_utilities.h

Build a local version:

python3 -m pip install -e .
# use 'CC="ccache clang" python3 -m pip install -e .' to speed up incremental builds

Publish

  1. Bump the version number in setup.py.

  2. Create a new release on GitHub.

MSVC compatibility (work in progress)