Merge pull request #2837 from bsipocz/irsa_sia

ENH: ipac.irsa to have SIA backend

Author: bsipocz

CHANGES.rst

@@ -26,6 +26,8 @@ ipac.irsa
   ``columns``, and the ``cache`` and ``verbose`` kwargs have been
   deprecated as they have no effect. [#2823]
 
+- Method to run SIAv2 VO queries is added. [#2837]
+
 gaia
 ^^^^
 
@@ -75,7 +77,7 @@ esa.xmm_newton
 simbad
 ^^^^^^
 
-- new ``query_tap`` method to access SIMBAD. This comes with additional methods to explore SIMBAD's tables and 
+- new ``query_tap`` method to access SIMBAD. This comes with additional methods to explore SIMBAD's tables and
   their links: ``Simbad.list_tables``, ``Simbad.list_columns``, and ``Simbad.list_linked_tables``. [#2856]
 
 solarsystem.neodys

astroquery/ipac/irsa/__init__.py

@@ -27,6 +27,7 @@ class Conf(_config.ConfigNamespace):
     timeout = _config.ConfigItem(
         60,
         'Time limit for connecting to the IRSA server.')
+    sia_url = _config.ConfigItem('https://irsa.ipac.caltech.edu/SIA', 'IRSA SIA URL')
     tap_url = _config.ConfigItem('https://irsa.ipac.caltech.edu/TAP', 'IRSA TAP URL')
 
 

astroquery/ipac/irsa/core.py

@@ -11,7 +11,16 @@
 from astropy.coordinates import SkyCoord, Angle
 from astropy import units as u
 from astropy.utils.decorators import deprecated_renamed_argument
+
 from pyvo.dal import TAPService
+
+try:
+    from pyvo.dal.sia2 import SIA2Service, SIA2_PARAMETERS_DESC
+except ImportError:
+    # Can be removed once min version of pyvo is 1.5
+    from pyvo.dal.sia2 import SIA_PARAMETERS_DESC as SIA2_PARAMETERS_DESC
+    from pyvo.dal.sia2 import SIAService as SIA2Service
+
 from astroquery import log
 from astroquery.query import BaseVOQuery
 from astroquery.utils.commons import parse_coordinates
@@ -26,9 +35,17 @@ class IrsaClass(BaseVOQuery):
 
     def __init__(self):
         super().__init__()
+        self.sia_url = conf.sia_url
         self.tap_url = conf.tap_url
+        self._sia = None
         self._tap = None
 
+    @property
+    def sia(self):
+        if not self._sia:
+            self._sia = SIA2Service(baseurl=self.sia_url, session=self._session)
+        return self._sia
+
     @property
     def tap(self):
         if not self._tap:
@@ -60,10 +77,59 @@ def query_tap(self, query, *, maxrec=None):
         log.debug(f'TAP query: {query}')
         return self.tap.search(query, language='ADQL', maxrec=maxrec)
 
+    def query_sia(self, *, pos=None, band=None, time=None, pol=None,
+                  field_of_view=None, spatial_resolution=None,
+                  spectral_resolving_power=None, exptime=None,
+                  timeres=None, publisher_did=None,
+                  facility=None, collection=None,
+                  instrument=None, data_type=None,
+                  calib_level=None, target_name=None,
+                  res_format=None, maxrec=None,
+                  **kwargs):
+        """
+        Use standard SIA2 attributes to query the IRSA SIA service.
+
+        Parameters
+        ----------
+        _SIA2_PARAMETERS
+
+        Returns
+        -------
+        Results in `pyvo.dal.SIAResults` format.
+        result.table in Astropy table format
+        """
+        return self.sia.search(
+            pos=pos,
+            band=band,
+            time=time,
+            pol=pol,
+            field_of_view=field_of_view,
+            spatial_resolution=spatial_resolution,
+            spectral_resolving_power=spectral_resolving_power,
+            exptime=exptime,
+            timeres=timeres,
+            publisher_did=publisher_did,
+            facility=facility,
+            collection=collection,
+            instrument=instrument,
+            data_type=data_type,
+            calib_level=calib_level,
+            target_name=target_name,
+            res_format=res_format,
+            maxrec=maxrec,
+            **kwargs)
+
+    # SIA2_PARAMETERS_DESC contains links that Sphinx can't resolve.
+    # SIA2_PARAMETERS_DESC contains links that Sphinx can't resolve.
+    for var in ('POLARIZATION_STATES', 'CALIBRATION_LEVELS'):
+        SIA2_PARAMETERS_DESC = SIA2_PARAMETERS_DESC.replace(f'`pyvo.dam.obscore.{var}`',
+                                                            f'pyvo.dam.obscore.{var}')
+    query_sia.__doc__ = query_sia.__doc__.replace('_SIA2_PARAMETERS', SIA2_PARAMETERS_DESC)
+
     @deprecated_renamed_argument(("selcols", "cache", "verbose"), ("columns", None, None), since="0.4.7")
     def query_region(self, coordinates=None, *, catalog=None, spatial='Cone',
                      radius=10 * u.arcsec, width=None, polygon=None,
-                     get_query_payload=False, columns=None,
+                     get_query_payload=False, columns='*',
                      verbose=False, cache=True):
         """
         Queries the IRSA TAP server around a coordinate and returns a `~astropy.table.Table` object.
@@ -106,9 +172,6 @@ def query_region(self, coordinates=None, *, catalog=None, spatial='Cone',
         if catalog is None:
             raise InvalidQueryError("Catalog name is required!")
 
-        if columns is None:
-            columns = '*'
-
         adql = f'SELECT {columns} FROM {catalog}'
 
         if spatial == 'All-Sky':

docs/ipac/irsa/irsa.rst

@@ -220,6 +220,63 @@ returned as a `~pyvo.dal.TAPResults` object. Its ``to_table`` or ``to_qtable`` m
     202.809023  46.964558  15.874     0.081 ...  15.322     0.188     AAC    000
 
 
+Simple image access queries
+---------------------------
+
+`~astroquery.ipac.irsa.IrsaClass.query_sia` provides a way to access IRSA's Simple
+Image Access VO service. In the following example we are looking for Spitzer
+Enhanced Imaging products in the centre of the COSMOS field as an `~astropy.table.Table`.
+
+.. doctest-remote-data::
+
+   >>> from astroquery.ipac.irsa import Irsa
+   >>> from astropy.coordinates import SkyCoord
+   >>> from astropy import units as u
+   >>>
+   >>> coord = SkyCoord('150.01d 2.2d', frame='icrs')
+   >>> spitzer_images = Irsa.query_sia(pos=(coord, 1 * u.arcmin), collection='spitzer_seip').to_table()
+
+Now open a cutout image for one of the science images. You could either use
+the the IRSA on-premise data or the cloud version of it using the
+``access_url`` or ``cloud_access`` columns. For more info about fits
+cutouts, please visit :ref:`astropy:fits_io_cloud`.
+
+.. doctest-remote-data::
+
+   >>> from astropy.io import fits
+   >>> from astropy.nddata import Cutout2D
+   >>> from astropy.wcs import WCS
+   >>> science_image = spitzer_images[spitzer_images['dataproduct_subtype'] == 'science'][0]
+   >>> with fits.open(science_image['access_url'], use_fsspec=True) as hdul:
+   ...     cutout = Cutout2D(hdul[0].section, position=coord, size=2 * u.arcmin, wcs=WCS(hdul[0].header))
+
+Now plot the cutout.
+
+.. doctest-skip::
+
+   >>> import matplotlib.pyplot as plt
+   >>> plt.imshow(cutout.data, cmap='grey')
+   >>> plt.show()
+
+.. plot::
+
+   from astroquery.ipac.irsa import Irsa
+   from astropy.coordinates import SkyCoord
+   from astropy import units as u
+   from astropy.io import fits
+   from astropy.nddata import Cutout2D
+   from astropy.wcs import WCS
+   import matplotlib.pyplot as plt
+   coord = SkyCoord('150.01d 2.2d', frame='icrs')
+   spitzer_images = Irsa.query_sia(pos=(coord, 1 * u.arcmin), collection='spitzer_seip').to_table()
+   science_image = spitzer_images[spitzer_images['dataproduct_subtype'] == 'science'][0]
+   with fits.open(science_image['access_url'], use_fsspec=True) as hdul:
+        cutout = Cutout2D(hdul[0].section, position=coord, size=2 *
+   u.arcmin, wcs=WCS(hdul[0].header))
+   plt.imshow(cutout.data, cmap='grey')
+   plt.show()
+
+
 Other Configurations
 --------------------
 

setup.cfg

@@ -152,6 +152,7 @@ docs=
    matplotlib
    sphinx-astropy>=1.5
    scipy
+   fsspec[http]
 all=
    mocpy>=0.9
    astropy-healpix