From d795cb7434b41c5dd959218151ebc8c71a1c6a7b Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Brigitta=20Sip=C5=91cz?= Date: Thu, 12 Sep 2024 17:52:29 -0700 Subject: [PATCH] Merge pull request #594 from duytnguyendtn/radquantity Support Astropy Quantity as radius arg for Registry SkyCoord Spatial constraint --- CHANGES.rst | 4 +--- docs/registry/index.rst | 31 +++++++++++++++++++++++------- pyvo/registry/rtcons.py | 15 +++++++++++++-- pyvo/registry/tests/test_rtcons.py | 9 +++++++++ 4 files changed, 47 insertions(+), 12 deletions(-) diff --git a/CHANGES.rst b/CHANGES.rst index 7363ddab..4385ed90 100644 --- a/CHANGES.rst +++ b/CHANGES.rst @@ -9,9 +9,7 @@ Bug Fixes - path separators are no longer taken over from image titles to file system paths. [#557] - -1.5.3 (unreleased) -================== +- Registry Spatial constraint now supports Astropy Quantities for the radius argument [#594] - Added `'sia1'` as servicetype for registry searches. [#583] diff --git a/docs/registry/index.rst b/docs/registry/index.rst index da1e3788..b8d1acfe 100644 --- a/docs/registry/index.rst +++ b/docs/registry/index.rst @@ -158,13 +158,30 @@ And to look for tap resources *in* a specific cone, you would do >>> from astropy.coordinates import SkyCoord >>> registry.search(registry.Freetext("Wolf-Rayet"), - ... registry.Spatial((SkyCoord("23d +3d"), 3), intersect="enclosed")) # doctest: +IGNORE_OUTPUT - - ivoid ... - ... - object ... - ---------------------------- ... - ivo://cds.vizier/j/aj/166/68 ... + ... registry.Spatial((SkyCoord("23d +3d"), 3), intersect="enclosed")) + + ivoid ... + ... + object ... + ------------------------------- ... + ivo://cds.vizier/j/a+a/688/a104 ... + ivo://cds.vizier/j/aj/166/68 ... + +Astropy Quantities are also supported for the radius angle of a SkyCoord-defined circular region: + +.. doctest-remote-data:: + + >>> from astropy.coordinates import SkyCoord + >>> from astropy import units as u + >>> registry.search(registry.Freetext("Wolf-Rayet"), + ... registry.Spatial((SkyCoord("23d +3d"), 180*u.Unit('arcmin')), intersect="enclosed")) + + ivoid ... + ... + object ... + ------------------------------- ... + ivo://cds.vizier/j/a+a/688/a104 ... + ivo://cds.vizier/j/aj/166/68 ... Where ``intersect`` can take the following values: * 'covers' is the default and returns resources that cover the geometry provided, diff --git a/pyvo/registry/rtcons.py b/pyvo/registry/rtcons.py index f41a763f..a07e6c9d 100644 --- a/pyvo/registry/rtcons.py +++ b/pyvo/registry/rtcons.py @@ -615,6 +615,10 @@ class Spatial(Constraint): are interpreted in degrees):: >>> resources = registry.Spatial((SkyCoord("23d +3d"), 3)) + + Or you can provide the radius angle as an Astropy Quantity: + + >>> resources = registry.Spatial((SkyCoord("23d +3d"), 1*u.rad)) """ _keyword = "spatial" _extra_tables = ["rr.stc_spatial"] @@ -631,7 +635,7 @@ def __init__(self, geom_spec, order=6, intersect="covers"): as a DALI point, a 3-sequence as a DALI circle, a 2n sequence as a DALI polygon. Additionally, strings are interpreted as ASCII MOCs, SkyCoords as points, and a pair of a - SkyCoord and a float as a circle. Other types (proper + SkyCoord and a float or Quantity as a circle. Other types (proper geometries or MOCPy objects) might be supported in the future. order : int, optional @@ -661,9 +665,16 @@ def tomoc(s): elif len(geom_spec) == 2: if isinstance(geom_spec[0], SkyCoord): + # If radius given is astropy quantity, then convert to degrees + if isinstance(geom_spec[1], u.Quantity): + if geom_spec[1].unit.physical_type != 'angle': + raise ValueError("Radius quantity is not of type angle.") + radius = geom_spec[1].to(u.deg).value + else: + radius = geom_spec[1] geom = tomoc(format_function_call("CIRCLE", [geom_spec[0].ra.value, geom_spec[0].dec.value, - geom_spec[1]])) + radius])) else: geom = tomoc(format_function_call("POINT", geom_spec)) diff --git a/pyvo/registry/tests/test_rtcons.py b/pyvo/registry/tests/test_rtcons.py index bf1f661f..60148ba2 100644 --- a/pyvo/registry/tests/test_rtcons.py +++ b/pyvo/registry/tests/test_rtcons.py @@ -253,6 +253,15 @@ def test_SkyCoord_Circle(self): assert cons.get_search_condition(FAKE_GAVO) == "1 = CONTAINS(MOC(6, CIRCLE(3.0, -30.0, 3)), coverage)" assert cons._extra_tables == ["rr.stc_spatial"] + def test_SkyCoord_Circle_RadiusQuantity(self): + for radius in [3*u.deg, 180*u.Unit('arcmin'), 10800*u.Unit('arcsec')]: + cons = registry.Spatial((SkyCoord(3 * u.deg, -30 * u.deg), radius)) + assert cons.get_search_condition(FAKE_GAVO) == ( + "1 = CONTAINS(MOC(6, CIRCLE(3.0, -30.0, 3.0)), coverage)") + + with pytest.raises(ValueError, match="is not of type angle."): + cons = registry.Spatial((SkyCoord(3 * u.deg, -30 * u.deg), (1 * u.m))) + def test_enclosed(self): cons = registry.Spatial("0/1-3", intersect="enclosed") assert cons.get_search_condition(FAKE_GAVO) == "1 = CONTAINS(coverage, MOC('0/1-3'))"