Skip to main content

target_match/
matcher.rs

1//! The matching engine: input trait, constraints, ranking, and a prebuilt index.
2//!
3//! A consumer's catalogue type implements [`SkyObject`] (position only — never a
4//! name). A [`Constraint`] combines a [`Membership`] shape (circular, rectangle,
5//! or rotated rectangle) with a [`Query`] mode. [`rank`] scans a slice; [`Matcher`]
6//! builds a declination-sorted index once and answers repeated queries, returning
7//! results identical to [`rank`].
8//!
9//! # Geometry
10//!
11//! Matching precesses the pointing to J2000 (via [`skymath::precess`]), then works
12//! in the local tangent frame about the pointing: an object's offset is decomposed
13//! from its great-circle separation and position angle (East of North, both from
14//! `skymath`) into East/North components, which are rotated into the camera frame
15//! for rectangular membership. The circumscribed circle pre-filters both rectangle
16//! tests, so the tangent decomposition is only evaluated for objects near the
17//! frame — never on the far side of the sky.
18
19use core::cmp::Ordering;
20
21use skymath::{position_angle, precess, separation, Angle, Epoch, Equatorial};
22
23use crate::optics::{Field, RadiusPolicy};
24
25/// A catalogue object that can be matched by sky position.
26///
27/// The trait exposes **only** a J2000 position — matching never reads a name or
28/// designation. A caller's own type keeps its identity; a [`Match`] borrows it.
29///
30/// # Example
31///
32/// ```
33/// use skymath::{Angle, Equatorial};
34/// use target_match::SkyObject;
35///
36/// struct Target {
37///     name: &'static str,
38///     ra: f64,
39///     dec: f64,
40/// }
41/// impl SkyObject for Target {
42///     fn position(&self) -> Equatorial {
43///         Equatorial::j2000(Angle::from_degrees(self.ra), Angle::from_degrees(self.dec)).unwrap()
44///     }
45/// }
46///
47/// let m31 = Target { name: "M 31", ra: 10.6847, dec: 41.2688 };
48/// assert!((m31.position().ra().degrees() - 10.6847).abs() < 1e-9);
49/// ```
50pub trait SkyObject {
51    /// The object's J2000 equatorial position.
52    fn position(&self) -> Equatorial;
53}
54
55/// The shape that decides whether an object is "in frame".
56///
57/// Usually built for you by a [`Constraint`] constructor; pass one directly to
58/// [`is_framed`] to test a single object.
59///
60/// # Example
61///
62/// ```
63/// use skymath::{Angle, Equatorial, ParseMode};
64/// use target_match::{is_framed, Membership, SkyObject};
65///
66/// struct Target {
67///     ra: f64,
68///     dec: f64,
69/// }
70/// impl SkyObject for Target {
71///     fn position(&self) -> Equatorial {
72///         Equatorial::j2000(Angle::from_degrees(self.ra), Angle::from_degrees(self.dec)).unwrap()
73///     }
74/// }
75///
76/// let m31 = Target { ra: 10.6847, dec: 41.2688 };
77/// let pointing = Equatorial::parse_j2000("00:42:44.3", "+41:16:09", ParseMode::Strict).unwrap();
78///
79/// let shape = Membership::Circular { radius: Angle::from_degrees(1.0) };
80/// assert!(is_framed(pointing, &m31, shape).in_frame);
81/// ```
82#[derive(Debug, Clone, Copy, PartialEq)]
83#[cfg_attr(feature = "serde", derive(serde::Serialize, serde::Deserialize))]
84pub enum Membership {
85    /// Pure angular distance: in frame iff separation ≤ `radius`.
86    Circular {
87        /// Search radius.
88        radius: Angle,
89    },
90    /// Axis-aligned rectangle of the given width×height field of view.
91    Rectangle {
92        /// `(width, height)` field of view.
93        fov: (Angle, Angle),
94    },
95    /// Rectangle rotated by a camera position angle (degrees East of North).
96    Rotated {
97        /// `(width, height)` field of view.
98        fov: (Angle, Angle),
99        /// Camera position angle, East of North.
100        position_angle: Angle,
101    },
102}
103
104/// What to return from a match.
105///
106/// Set on a [`Constraint`] via [`Constraint::all`], [`Constraint::nearest_one`],
107/// [`Constraint::nearest_n`], or [`Constraint::nearest_n_within`].
108///
109/// # Example
110///
111/// ```
112/// use skymath::{Angle, Equatorial, ParseMode};
113/// use target_match::{rank, Constraint, SkyObject};
114///
115/// struct Target {
116///     ra: f64,
117///     dec: f64,
118/// }
119/// impl SkyObject for Target {
120///     fn position(&self) -> Equatorial {
121///         Equatorial::j2000(Angle::from_degrees(self.ra), Angle::from_degrees(self.dec)).unwrap()
122///     }
123/// }
124///
125/// let catalog = [Target { ra: 10.6847, dec: 41.2688 }];
126/// let pointing = Equatorial::parse_j2000("00:42:44.3", "+41:16:09", ParseMode::Strict).unwrap();
127///
128/// // `nearest_n(1)` sets the query to `Query::NearestN { n: 1, .. }`.
129/// let hits = rank(pointing, &catalog, Constraint::circular(Angle::from_degrees(1.0)).nearest_n(1));
130/// assert_eq!(hits.len(), 1);
131/// ```
132#[derive(Debug, Clone, Copy, PartialEq)]
133#[cfg_attr(feature = "serde", derive(serde::Serialize, serde::Deserialize))]
134pub enum Query {
135    /// Every object inside the membership shape, ranked nearest-first.
136    AllWithinField,
137    /// The single nearest in-frame object.
138    NearestOne,
139    /// The `n` nearest objects by separation, optionally bounded by a radius.
140    NearestN {
141        /// Maximum number of results.
142        n: usize,
143        /// Optional maximum separation; unbounded when `None`.
144        max_radius: Option<Angle>,
145    },
146}
147
148/// A [`Membership`] shape combined with a [`Query`] mode (and the plate scale,
149/// when known, so pixel offsets can be reported).
150///
151/// Build one with a shape constructor ([`within`](Constraint::within),
152/// [`circular`](Constraint::circular), [`frame`](Constraint::frame),
153/// [`frame_rotated`](Constraint::frame_rotated)), then optionally switch the
154/// query mode ([`all`](Constraint::all), [`nearest_one`](Constraint::nearest_one),
155/// [`nearest_n`](Constraint::nearest_n),
156/// [`nearest_n_within`](Constraint::nearest_n_within)). Pass the result to
157/// [`rank`], [`Matcher::query`], or [`is_framed`].
158///
159/// # Example
160///
161/// ```
162/// use skymath::{Angle, Equatorial, ParseMode};
163/// use target_match::{rank, Constraint, Field, Optics, RadiusPolicy, SkyObject};
164///
165/// struct Target {
166///     ra: f64,
167///     dec: f64,
168/// }
169/// impl SkyObject for Target {
170///     fn position(&self) -> Equatorial {
171///         Equatorial::j2000(Angle::from_degrees(self.ra), Angle::from_degrees(self.dec)).unwrap()
172///     }
173/// }
174///
175/// let catalog = [Target { ra: 10.6847, dec: 41.2688 }];
176/// let pointing = Equatorial::parse_j2000("00:42:44.3", "+41:16:09", ParseMode::Strict).unwrap();
177/// let field = Field::from_optics(Optics {
178///     focal_mm: 800.0,
179///     pixel_um: (3.76, 3.76),
180///     binning: (1, 1),
181///     pixels: (6248, 4176),
182/// })
183/// .unwrap();
184///
185/// let c = Constraint::within(&field, RadiusPolicy::Circumscribed).nearest_one();
186/// assert_eq!(rank(pointing, &catalog, c).len(), 1);
187/// ```
188#[derive(Debug, Clone, Copy, PartialEq)]
189#[cfg_attr(feature = "serde", derive(serde::Serialize, serde::Deserialize))]
190pub struct Constraint {
191    /// The in-frame shape.
192    pub membership: Membership,
193    /// The query mode.
194    pub query: Query,
195    /// Per-axis plate scale (arcsec/px), used only to fill pixel offsets.
196    pub pixel_scale: Option<(f64, f64)>,
197}
198
199impl Constraint {
200    /// Circular membership with a radius derived from a field under `policy`.
201    ///
202    /// # Example
203    ///
204    /// ```
205    /// use target_match::{Constraint, Field, Optics, RadiusPolicy};
206    ///
207    /// let field = Field::from_optics(Optics {
208    ///     focal_mm: 800.0,
209    ///     pixel_um: (3.76, 3.76),
210    ///     binning: (1, 1),
211    ///     pixels: (6248, 4176),
212    /// })
213    /// .unwrap();
214    /// let c = Constraint::within(&field, RadiusPolicy::Circumscribed);
215    /// assert!(c.pixel_scale.is_some(), "the field's plate scale carries through");
216    /// ```
217    #[must_use]
218    pub fn within(field: &Field, policy: RadiusPolicy) -> Self {
219        Self {
220            membership: Membership::Circular {
221                radius: field.radius(policy),
222            },
223            query: Query::AllWithinField,
224            pixel_scale: field.pixel_scale(),
225        }
226    }
227    /// Circular membership with an explicit radius (no plate scale).
228    ///
229    /// # Example
230    ///
231    /// ```
232    /// use skymath::Angle;
233    /// use target_match::Constraint;
234    ///
235    /// let c = Constraint::circular(Angle::from_degrees(2.0));
236    /// assert!(c.pixel_scale.is_none(), "no `Field`, so no plate scale");
237    /// ```
238    #[must_use]
239    pub fn circular(radius: Angle) -> Self {
240        Self {
241            membership: Membership::Circular { radius },
242            query: Query::AllWithinField,
243            pixel_scale: None,
244        }
245    }
246    /// Axis-aligned rectangular membership from a field's width×height.
247    ///
248    /// # Example
249    ///
250    /// ```
251    /// use target_match::{Constraint, Field, Membership, Optics};
252    ///
253    /// let field = Field::from_optics(Optics {
254    ///     focal_mm: 800.0,
255    ///     pixel_um: (3.76, 3.76),
256    ///     binning: (1, 1),
257    ///     pixels: (6248, 4176),
258    /// })
259    /// .unwrap();
260    /// let c = Constraint::frame(&field);
261    /// assert!(matches!(c.membership, Membership::Rectangle { .. }));
262    /// ```
263    #[must_use]
264    pub fn frame(field: &Field) -> Self {
265        Self {
266            membership: Membership::Rectangle {
267                fov: (field.width(), field.height()),
268            },
269            query: Query::AllWithinField,
270            pixel_scale: field.pixel_scale(),
271        }
272    }
273    /// Rotated rectangular membership from a field and a camera position angle
274    /// (degrees East of North).
275    ///
276    /// # Example
277    ///
278    /// An object due north of the pointing sits on the frame's `+y` axis when
279    /// axis-aligned; a 90° East-of-North camera rotation moves it onto the `-x`
280    /// axis.
281    ///
282    /// ```
283    /// use skymath::{Angle, Equatorial, ParseMode};
284    /// use target_match::{rank, Constraint, Field, SkyObject};
285    ///
286    /// struct Target {
287    ///     ra: f64,
288    ///     dec: f64,
289    /// }
290    /// impl SkyObject for Target {
291    ///     fn position(&self) -> Equatorial {
292    ///         Equatorial::j2000(Angle::from_degrees(self.ra), Angle::from_degrees(self.dec)).unwrap()
293    ///     }
294    /// }
295    ///
296    /// let pointing = Equatorial::parse_j2000("00:42:44.3", "+41:16:09", ParseMode::Strict).unwrap();
297    /// let field = Field::from_fov(Angle::from_degrees(2.0), Angle::from_degrees(2.0)).unwrap();
298    /// let catalog = [Target { ra: 10.6847, dec: 41.2688 + 0.3 }]; // ~0.3° due north
299    ///
300    /// let (ax, ay) = rank(pointing, &catalog, Constraint::frame(&field))[0].offset.frame.unwrap();
301    /// assert!(ay.degrees() > 0.25 && ax.degrees().abs() < 0.05);
302    ///
303    /// let rotated = Constraint::frame_rotated(&field, Angle::from_degrees(90.0));
304    /// let (rx, ry) = rank(pointing, &catalog, rotated)[0].offset.frame.unwrap();
305    /// assert!(rx.degrees() < -0.25 && ry.degrees().abs() < 0.05);
306    /// ```
307    #[must_use]
308    pub fn frame_rotated(field: &Field, position_angle: Angle) -> Self {
309        Self {
310            membership: Membership::Rotated {
311                fov: (field.width(), field.height()),
312                position_angle,
313            },
314            query: Query::AllWithinField,
315            pixel_scale: field.pixel_scale(),
316        }
317    }
318    /// Set the query to all-within-field.
319    ///
320    /// # Example
321    ///
322    /// ```
323    /// use skymath::Angle;
324    /// use target_match::{Constraint, Query};
325    ///
326    /// let c = Constraint::circular(Angle::from_degrees(1.0)).nearest_one().all();
327    /// assert_eq!(c.query, Query::AllWithinField);
328    /// ```
329    #[must_use]
330    pub fn all(mut self) -> Self {
331        self.query = Query::AllWithinField;
332        self
333    }
334    /// Set the query to nearest-one.
335    ///
336    /// # Example
337    ///
338    /// ```
339    /// use skymath::Angle;
340    /// use target_match::{Constraint, Query};
341    ///
342    /// let c = Constraint::circular(Angle::from_degrees(1.0)).nearest_one();
343    /// assert_eq!(c.query, Query::NearestOne);
344    /// ```
345    #[must_use]
346    pub fn nearest_one(mut self) -> Self {
347        self.query = Query::NearestOne;
348        self
349    }
350    /// Set the query to the `n` nearest (unbounded).
351    ///
352    /// # Example
353    ///
354    /// ```
355    /// use skymath::Angle;
356    /// use target_match::{Constraint, Query};
357    ///
358    /// let c = Constraint::circular(Angle::from_degrees(1.0)).nearest_n(3);
359    /// assert_eq!(c.query, Query::NearestN { n: 3, max_radius: None });
360    /// ```
361    #[must_use]
362    pub fn nearest_n(mut self, n: usize) -> Self {
363        self.query = Query::NearestN {
364            n,
365            max_radius: None,
366        };
367        self
368    }
369    /// Set the query to the `n` nearest within `max_radius`.
370    ///
371    /// # Example
372    ///
373    /// ```
374    /// use skymath::Angle;
375    /// use target_match::{Constraint, Query};
376    ///
377    /// let radius = Angle::from_degrees(1.0);
378    /// let c = Constraint::circular(radius).nearest_n_within(3, radius);
379    /// assert_eq!(c.query, Query::NearestN { n: 3, max_radius: Some(radius) });
380    /// ```
381    #[must_use]
382    pub fn nearest_n_within(mut self, n: usize, max_radius: Angle) -> Self {
383        self.query = Query::NearestN {
384            n,
385            max_radius: Some(max_radius),
386        };
387        self
388    }
389}
390
391/// The offset of a matched object relative to the frame centre.
392///
393/// Carried on every [`Match`]. `frame` and `pixels` are only populated for
394/// rectangular [`Membership`] (see [`Constraint::frame`],
395/// [`Constraint::frame_rotated`]) — circular membership has no frame axes.
396///
397/// # Example
398///
399/// ```
400/// use skymath::{Angle, Equatorial, ParseMode};
401/// use target_match::{rank, Constraint, Field, Optics, SkyObject};
402///
403/// struct Target {
404///     ra: f64,
405///     dec: f64,
406/// }
407/// impl SkyObject for Target {
408///     fn position(&self) -> Equatorial {
409///         Equatorial::j2000(Angle::from_degrees(self.ra), Angle::from_degrees(self.dec)).unwrap()
410///     }
411/// }
412///
413/// let catalog = [Target { ra: 10.6847, dec: 41.2688 }];
414/// let pointing = Equatorial::parse_j2000("00:42:44.3", "+41:16:09", ParseMode::Strict).unwrap();
415/// let field = Field::from_optics(Optics {
416///     focal_mm: 800.0,
417///     pixel_um: (3.76, 3.76),
418///     binning: (1, 1),
419///     pixels: (6248, 4176),
420/// })
421/// .unwrap();
422///
423/// let hits = rank(pointing, &catalog, Constraint::frame(&field));
424/// let offset = hits[0].offset;
425/// assert!(offset.frame.is_some(), "rectangular membership reports a frame-aligned offset");
426/// assert!(offset.pixels.is_some(), "plate scale from `Optics` fills the pixel offset");
427/// ```
428#[derive(Debug, Clone, Copy, PartialEq)]
429#[cfg_attr(feature = "serde", derive(serde::Serialize, serde::Deserialize))]
430pub struct Offset {
431    /// Sky-tangent offset `(East, North)` — always present.
432    pub sky: (Angle, Angle),
433    /// Frame-aligned offset `(x, y)`, present for rectangular membership.
434    pub frame: Option<(Angle, Angle)>,
435    /// Frame-aligned offset in pixels `(x, y)`, present when a plate scale is known.
436    pub pixels: Option<(f64, f64)>,
437}
438
439/// A ranked match: a borrowed catalogue object plus its computed geometry.
440///
441/// Returned by [`rank`], [`Matcher::query`], and [`is_framed`].
442///
443/// # Example
444///
445/// ```
446/// use skymath::{Angle, Equatorial, ParseMode};
447/// use target_match::{rank, Constraint, Field, Optics, RadiusPolicy, SkyObject};
448///
449/// struct Target {
450///     name: &'static str,
451///     ra: f64,
452///     dec: f64,
453/// }
454/// impl SkyObject for Target {
455///     fn position(&self) -> Equatorial {
456///         Equatorial::j2000(Angle::from_degrees(self.ra), Angle::from_degrees(self.dec)).unwrap()
457///     }
458/// }
459///
460/// let catalog = [Target { name: "M 31", ra: 10.6847, dec: 41.2688 }];
461/// let pointing = Equatorial::parse_j2000("00:42:44.3", "+41:16:09", ParseMode::Strict).unwrap();
462/// let field = Field::from_optics(Optics {
463///     focal_mm: 800.0,
464///     pixel_um: (3.76, 3.76),
465///     binning: (1, 1),
466///     pixels: (6248, 4176),
467/// })
468/// .unwrap();
469///
470/// let hits = rank(pointing, &catalog, Constraint::within(&field, RadiusPolicy::Circumscribed));
471/// let m = &hits[0];
472/// assert_eq!(m.object.name, "M 31");
473/// assert!(m.in_frame);
474/// assert!(m.separation.arcseconds() < 2.0);
475/// ```
476#[derive(Debug, Clone, Copy, PartialEq)]
477pub struct Match<'a, T> {
478    /// The matched object, borrowed from the caller's slice or the [`Matcher`].
479    pub object: &'a T,
480    /// Great-circle separation from the frame centre.
481    pub separation: Angle,
482    /// Whether the object is inside the active membership shape.
483    pub in_frame: bool,
484    /// The object's offset relative to the frame centre.
485    pub offset: Offset,
486    /// Position angle from frame centre to the object (degrees East of North).
487    pub position_angle: Angle,
488}
489
490/// Tangent-frame `(East, North)` offset in radians, decomposed from an
491/// already-computed separation and position angle (the polar decomposition
492/// `skymath::tangent_offset` performs, reusing this evaluation's sep/PA
493/// instead of recomputing them).
494fn tangent_components(sep: Angle, pa: Angle) -> (f64, f64) {
495    let (s, r) = (pa.radians(), sep.radians());
496    (r * s.sin(), r * s.cos())
497}
498
499fn rotate(east: f64, north: f64, pa_rad: f64) -> (f64, f64) {
500    let (s, c) = pa_rad.sin_cos();
501    (east * c - north * s, east * s + north * c)
502}
503
504/// The circumscribed-circle radius (radians) that bounds a membership shape.
505fn bound_radius(m: Membership) -> f64 {
506    match m {
507        Membership::Circular { radius } => radius.radians(),
508        Membership::Rectangle { fov } | Membership::Rotated { fov, .. } => {
509            (fov.0.radians() / 2.0).hypot(fov.1.radians() / 2.0)
510        }
511    }
512}
513
514fn contains(sep: Angle, east: f64, north: f64, m: Membership) -> bool {
515    match m {
516        Membership::Circular { radius } => {
517            let r = radius.radians();
518            r.is_finite() && r >= 0.0 && sep.radians() <= r
519        }
520        Membership::Rectangle { fov } => within_rect(sep, east, north, fov, 0.0),
521        Membership::Rotated {
522            fov,
523            position_angle,
524        } => within_rect(sep, east, north, fov, position_angle.radians()),
525    }
526}
527
528fn within_rect(sep: Angle, east: f64, north: f64, fov: (Angle, Angle), pa: f64) -> bool {
529    let (hx, hy) = (fov.0.radians() / 2.0, fov.1.radians() / 2.0);
530    let circum = hx.hypot(hy);
531    // Circumscribed pre-filter (also rejects NaN separations).
532    if sep.radians() > circum || sep.radians().is_nan() {
533        return false;
534    }
535    let (x, y) = rotate(east, north, pa);
536    x.abs() <= hx && y.abs() <= hy
537}
538
539fn build_offset(east: f64, north: f64, m: Membership, scale: Option<(f64, f64)>) -> Offset {
540    let sky = (Angle::from_radians(east), Angle::from_radians(north));
541    match m {
542        Membership::Circular { .. } => Offset {
543            sky,
544            frame: None,
545            pixels: None,
546        },
547        Membership::Rectangle { .. } | Membership::Rotated { .. } => {
548            let pa = match m {
549                Membership::Rotated { position_angle, .. } => position_angle.radians(),
550                _ => 0.0,
551            };
552            let (x, y) = rotate(east, north, pa);
553            let frame = Some((Angle::from_radians(x), Angle::from_radians(y)));
554            let pixels = scale.map(|(sx, sy)| {
555                (
556                    Angle::from_radians(x).arcseconds() / sx,
557                    Angle::from_radians(y).arcseconds() / sy,
558                )
559            });
560            Offset { sky, frame, pixels }
561        }
562    }
563}
564
565fn evaluate<'a, T: SkyObject>(
566    pointing: Equatorial,
567    obj: &'a T,
568    m: Membership,
569    scale: Option<(f64, f64)>,
570) -> Match<'a, T> {
571    let pos = obj.position();
572    let sep = separation(pointing, pos);
573    let pa = position_angle(pointing, pos);
574    let (east, north) = tangent_components(sep, pa);
575    Match {
576        object: obj,
577        separation: sep,
578        in_frame: contains(sep, east, north, m),
579        offset: build_offset(east, north, m, scale),
580        position_angle: pa,
581    }
582}
583
584/// Whether `obj` should be kept for `query` given its evaluated match.
585fn keep<T>(m: &Match<'_, T>, query: Query) -> bool {
586    match query {
587        Query::AllWithinField | Query::NearestOne => m.in_frame,
588        Query::NearestN { max_radius, .. } => {
589            max_radius.map_or(true, |r| m.separation.radians() <= r.radians())
590        }
591    }
592}
593
594/// Shared core: evaluate candidates, filter, rank, and truncate per the query.
595fn rank_candidates<'a, T: SkyObject, I>(
596    pointing: Equatorial,
597    candidates: I,
598    c: &Constraint,
599) -> Vec<Match<'a, T>>
600where
601    I: Iterator<Item = (usize, &'a T)>,
602{
603    let mut scored: Vec<(usize, Match<'a, T>)> = candidates
604        .map(|(i, o)| (i, evaluate(pointing, o, c.membership, c.pixel_scale)))
605        .filter(|(_, m)| keep(m, c.query))
606        .collect();
607    scored.sort_by(|a, b| {
608        a.1.separation
609            .radians()
610            .partial_cmp(&b.1.separation.radians())
611            .unwrap_or(Ordering::Equal)
612            .then(a.0.cmp(&b.0))
613    });
614    let mut out: Vec<Match<'a, T>> = scored.into_iter().map(|(_, m)| m).collect();
615    match c.query {
616        Query::NearestOne => out.truncate(1),
617        Query::NearestN { n, .. } => out.truncate(n),
618        Query::AllWithinField => {}
619    }
620    out
621}
622
623/// Rank a slice of objects against a pointing under a constraint (stateless scan).
624///
625/// The pointing is precessed to J2000 first. Results are ascending by separation
626/// with ties broken by input order. For repeated queries against one catalogue,
627/// build a [`Matcher`] instead — it returns identical results faster.
628///
629/// # Example
630///
631/// ```
632/// use skymath::{Angle, Equatorial, ParseMode};
633/// use target_match::{rank, Constraint, Field, Optics, RadiusPolicy, SkyObject};
634///
635/// struct Target {
636///     name: &'static str,
637///     ra: f64,
638///     dec: f64,
639/// }
640/// impl SkyObject for Target {
641///     fn position(&self) -> Equatorial {
642///         Equatorial::j2000(Angle::from_degrees(self.ra), Angle::from_degrees(self.dec)).unwrap()
643///     }
644/// }
645///
646/// let catalog = [
647///     Target { name: "M 31", ra: 10.6847, dec: 41.2688 },
648///     Target { name: "M 33", ra: 23.4621, dec: 30.6599 },
649/// ];
650/// let pointing = Equatorial::parse_j2000("00:42:44.3", "+41:16:09", ParseMode::Strict).unwrap();
651/// let field = Field::from_optics(Optics {
652///     focal_mm: 800.0,
653///     pixel_um: (3.76, 3.76),
654///     binning: (1, 1),
655///     pixels: (6248, 4176),
656/// })
657/// .unwrap();
658///
659/// let hits = rank(pointing, &catalog, Constraint::within(&field, RadiusPolicy::Circumscribed).nearest_one());
660/// assert_eq!(hits[0].object.name, "M 31");
661/// ```
662#[must_use]
663pub fn rank<T: SkyObject>(pointing: Equatorial, objects: &[T], c: Constraint) -> Vec<Match<'_, T>> {
664    let p = precess(pointing, Epoch::J2000);
665    rank_candidates(p, objects.iter().enumerate(), &c)
666}
667
668/// Evaluate a single object against a frame, returning its membership + geometry.
669///
670/// The pointing is precessed to J2000 first. Unlike [`rank`], no filtering or
671/// ranking is applied — the returned [`Match`] always describes `object`.
672///
673/// # Example
674///
675/// ```
676/// use skymath::{Angle, Equatorial, ParseMode};
677/// use target_match::{is_framed, Membership, SkyObject};
678///
679/// struct Target {
680///     ra: f64,
681///     dec: f64,
682/// }
683/// impl SkyObject for Target {
684///     fn position(&self) -> Equatorial {
685///         Equatorial::j2000(Angle::from_degrees(self.ra), Angle::from_degrees(self.dec)).unwrap()
686///     }
687/// }
688///
689/// let m31 = Target { ra: 10.6847, dec: 41.2688 };
690/// let pointing = Equatorial::parse_j2000("00:42:44.3", "+41:16:09", ParseMode::Strict).unwrap();
691///
692/// let m = is_framed(pointing, &m31, Membership::Circular { radius: Angle::from_degrees(1.0) });
693/// assert!(m.in_frame);
694/// ```
695#[must_use]
696pub fn is_framed<T: SkyObject>(
697    pointing: Equatorial,
698    object: &T,
699    membership: Membership,
700) -> Match<'_, T> {
701    let p = precess(pointing, Epoch::J2000);
702    evaluate(p, object, membership, None)
703}
704
705/// A prebuilt, declination-sorted index for repeated queries against one catalogue.
706///
707/// Produces results identical to [`rank`] for the same objects, pointing, and
708/// constraint — the index is a performance optimization only. Build it once, then
709/// [`query`](Matcher::query) many pointings.
710///
711/// # Example
712///
713/// ```
714/// use skymath::{Angle, Equatorial, ParseMode};
715/// use target_match::{Constraint, Field, Matcher, Membership, Optics, RadiusPolicy, SkyObject};
716///
717/// struct Target {
718///     name: &'static str,
719///     ra: f64,
720///     dec: f64,
721/// }
722/// impl SkyObject for Target {
723///     fn position(&self) -> Equatorial {
724///         Equatorial::j2000(Angle::from_degrees(self.ra), Angle::from_degrees(self.dec)).unwrap()
725///     }
726/// }
727///
728/// let matcher = Matcher::from_objects(vec![
729///     Target { name: "M 31", ra: 10.6847, dec: 41.2688 },
730///     Target { name: "M 33", ra: 23.4621, dec: 30.6599 },
731/// ]);
732/// assert_eq!(matcher.objects().len(), 2);
733///
734/// let pointing = Equatorial::parse_j2000("00:42:44.3", "+41:16:09", ParseMode::Strict).unwrap();
735/// let field = Field::from_optics(Optics {
736///     focal_mm: 800.0,
737///     pixel_um: (3.76, 3.76),
738///     binning: (1, 1),
739///     pixels: (6248, 4176),
740/// })
741/// .unwrap();
742///
743/// let hits = matcher.query(pointing, Constraint::within(&field, RadiusPolicy::Circumscribed).nearest_one());
744/// assert_eq!(hits[0].object.name, "M 31");
745///
746/// let m = matcher.is_framed(
747///     pointing,
748///     &matcher.objects()[0],
749///     Membership::Circular { radius: Angle::from_degrees(1.0) },
750/// );
751/// assert!(m.in_frame);
752/// ```
753pub struct Matcher<T> {
754    storage: Vec<T>,
755    /// `(declination_degrees, original_index)` sorted by declination.
756    sorted: Vec<(f64, usize)>,
757}
758
759impl<T: SkyObject> Matcher<T> {
760    /// Build an index from a set of objects (original order is preserved for
761    /// tie-breaking and [`objects`](Matcher::objects)).
762    #[must_use]
763    pub fn from_objects(objects: Vec<T>) -> Self {
764        let mut sorted: Vec<(f64, usize)> = objects
765            .iter()
766            .enumerate()
767            .map(|(i, o)| (o.position().dec().degrees(), i))
768            .collect();
769        sorted.sort_by(|a, b| {
770            a.0.partial_cmp(&b.0)
771                .unwrap_or(Ordering::Equal)
772                .then(a.1.cmp(&b.1))
773        });
774        Self {
775            storage: objects,
776            sorted,
777        }
778    }
779
780    /// The stored objects, in their original insertion order.
781    #[must_use]
782    pub fn objects(&self) -> &[T] {
783        &self.storage
784    }
785
786    /// Query the index for a pointing under a constraint.
787    #[must_use]
788    pub fn query(&self, pointing: Equatorial, c: Constraint) -> Vec<Match<'_, T>> {
789        let p = precess(pointing, Epoch::J2000);
790        let r = match c.query {
791            Query::NearestN { max_radius, .. } => max_radius.map_or(f64::INFINITY, |a| a.radians()),
792            _ => bound_radius(c.membership),
793        };
794        let idxs = self.band(p.dec().degrees(), r);
795        rank_candidates(p, idxs.into_iter().map(|i| (i, &self.storage[i])), &c)
796    }
797
798    /// Evaluate a single stored-or-external object against a frame (see [`is_framed`]).
799    #[must_use]
800    pub fn is_framed<'a>(
801        &self,
802        pointing: Equatorial,
803        object: &'a T,
804        m: Membership,
805    ) -> Match<'a, T> {
806        is_framed(pointing, object, m)
807    }
808
809    /// Original indices whose declination lies within `r_rad` of `dec0_deg`.
810    fn band(&self, dec0_deg: f64, r_rad: f64) -> Vec<usize> {
811        if r_rad.is_infinite() && r_rad > 0.0 {
812            return self.sorted.iter().map(|&(_, i)| i).collect();
813        }
814        if !r_rad.is_finite() || r_rad < 0.0 {
815            return Vec::new();
816        }
817        let r_deg = r_rad.to_degrees();
818        let (lo, hi) = (dec0_deg - r_deg, dec0_deg + r_deg);
819        let start = self.sorted.partition_point(|&(d, _)| d < lo);
820        let end = self.sorted.partition_point(|&(d, _)| d <= hi);
821        self.sorted[start..end].iter().map(|&(_, i)| i).collect()
822    }
823}
824
825#[cfg(test)]
826mod tests {
827    use super::*;
828
829    #[derive(Clone)]
830    struct Obj {
831        name: &'static str,
832        ra: f64,
833        dec: f64,
834    }
835    impl SkyObject for Obj {
836        fn position(&self) -> Equatorial {
837            Equatorial::j2000(Angle::from_degrees(self.ra), Angle::from_degrees(self.dec)).unwrap()
838        }
839    }
840
841    fn catalog() -> Vec<Obj> {
842        vec![
843            Obj {
844                name: "M 31",
845                ra: 10.6847,
846                dec: 41.2688,
847            },
848            Obj {
849                name: "M 110",
850                ra: 10.0921,
851                dec: 41.6853,
852            },
853            Obj {
854                name: "M 33",
855                ra: 23.4621,
856                dec: 30.6599,
857            },
858            Obj {
859                name: "M 42",
860                ra: 83.8221,
861                dec: -5.3911,
862            },
863        ]
864    }
865
866    fn m31() -> Equatorial {
867        Equatorial::j2000(Angle::from_degrees(10.6847), Angle::from_degrees(41.2688)).unwrap()
868    }
869
870    #[test]
871    fn nearest_one_circular() {
872        let cat = catalog();
873        let c = Constraint::circular(Angle::from_degrees(2.0)).nearest_one();
874        let hits = rank(m31(), &cat, c);
875        assert_eq!(hits.len(), 1);
876        assert_eq!(hits[0].object.name, "M 31");
877        assert!(hits[0].separation.arcseconds() < 1.0);
878    }
879
880    #[test]
881    fn all_within_field_ranked() {
882        let cat = catalog();
883        // ~1° radius keeps M31 (self) + M110 (~0.62°); excludes M33/M42.
884        let c = Constraint::circular(Angle::from_degrees(1.0)).all();
885        let hits = rank(m31(), &cat, c);
886        assert_eq!(
887            hits.len(),
888            2,
889            "{hits:?}",
890            hits = hits.iter().map(|h| h.object.name).collect::<Vec<_>>()
891        );
892        assert_eq!(hits[0].object.name, "M 31");
893        assert_eq!(hits[1].object.name, "M 110");
894        assert!(hits[0].separation.radians() <= hits[1].separation.radians());
895    }
896
897    #[test]
898    fn nearest_n_bounds_and_counts() {
899        let cat = catalog();
900        let c = Constraint::circular(Angle::from_degrees(1.0)).nearest_n(3);
901        let hits = rank(m31(), &cat, c);
902        assert_eq!(hits.len(), 3, "top-3 by separation regardless of frame");
903        assert_eq!(hits[0].object.name, "M 31");
904        // Bounded nearest-N respects the radius.
905        let c2 = Constraint::circular(Angle::from_degrees(1.0))
906            .nearest_n_within(3, Angle::from_degrees(1.0));
907        assert_eq!(rank(m31(), &cat, c2).len(), 2, "only M31 + M110 within 1°");
908    }
909
910    #[test]
911    fn coordinates_only_never_name() {
912        // A far object literally named "M 31" must NOT match near M31's pointing.
913        let cat = vec![
914            Obj {
915                name: "M 31",
916                ra: 200.0,
917                dec: -40.0,
918            },
919            Obj {
920                name: "Some Galaxy",
921                ra: 10.6847,
922                dec: 41.2688,
923            },
924        ];
925        let c = Constraint::circular(Angle::from_degrees(2.0)).nearest_one();
926        let hits = rank(m31(), &cat, c);
927        assert_eq!(hits.len(), 1);
928        assert_eq!(hits[0].object.name, "Some Galaxy");
929    }
930
931    #[test]
932    fn rectangle_excludes_circle_only_corner() {
933        // An object just outside the rectangle but inside the circumscribed circle.
934        // Field 2°×1°: half-width 1°, half-height 0.5°. Put an object 0.9° North
935        // (in frame) vs 0.9° along the diagonal (out of the axis-aligned rect).
936        let field = Field::from_fov(Angle::from_degrees(2.0), Angle::from_degrees(1.0)).unwrap();
937        let north_obj = Obj {
938            name: "N",
939            ra: 10.6847,
940            dec: 41.2688 + 0.4,
941        }; // within 0.5° height
942        let high_obj = Obj {
943            name: "H",
944            ra: 10.6847,
945            dec: 41.2688 + 0.9,
946        }; // beyond height, inside circle
947        let cat = vec![north_obj, high_obj];
948        let c = Constraint::frame(&field).all();
949        let hits = rank(m31(), &cat, c);
950        assert_eq!(hits.len(), 1);
951        assert_eq!(hits[0].object.name, "N");
952    }
953
954    #[test]
955    fn matcher_matches_rank_exactly() {
956        let cat = catalog();
957        let c = Constraint::circular(Angle::from_degrees(5.0)).all();
958        let via_rank: Vec<_> = rank(m31(), &cat, c).iter().map(|m| m.object.name).collect();
959        let matcher = Matcher::from_objects(cat.clone());
960        let via_index: Vec<_> = matcher
961            .query(m31(), c)
962            .iter()
963            .map(|m| m.object.name)
964            .collect();
965        assert_eq!(via_rank, via_index);
966        assert_eq!(matcher.objects().len(), 4);
967    }
968
969    #[test]
970    fn is_framed_reports_geometry() {
971        let m110 = catalog()[1].clone();
972        let m = is_framed(
973            m31(),
974            &m110,
975            Membership::Circular {
976                radius: Angle::from_degrees(1.0),
977            },
978        );
979        assert!(m.in_frame);
980        assert!((0.4..0.9).contains(&m.separation.degrees()));
981    }
982
983    #[test]
984    fn empty_catalog_and_zero_radius() {
985        let cat = catalog();
986        assert!(rank(
987            m31(),
988            &[] as &[Obj],
989            Constraint::circular(Angle::from_degrees(1.0))
990        )
991        .is_empty());
992        let c = Constraint::circular(Angle::from_degrees(-1.0)).all();
993        assert!(
994            rank(m31(), &cat, c).is_empty(),
995            "negative radius matches nothing"
996        );
997    }
998}