Desi legacy survey

astroquery/desi/__init__.py

@@ -0,0 +1,42 @@
+# Licensed under a 3-clause BSD style license - see LICENSE.rst
+
+"""
+DESI LegacySurvery
+
+https://www.legacysurvey.org/
+-----------------------------
+
+:author: Gabriele Barni (Gabriele.Barni@unige.ch)
+"""
+
+# Make the URL of the server, timeout and other items configurable
+# See <http://docs.astropy.org/en/latest/config/index.html#developer-usage>
+# for docs and examples on how to do this
+# Below is a common use case
+from astropy import config as _config
+
+
+class Conf(_config.ConfigNamespace):
+    """
+    Configuration parameters for `astroquery.desi`.
+    """
+    legacysurvey_service_url = _config.ConfigItem(
+        ['https://www.legacysurvey.org/viewer/fits-cutout',
+         ],
+        'url for the LegacySurvey service')
+
+    tap_service_url = _config.ConfigItem(
+        ['https://datalab.noirlab.edu/tap',
+         ],
+        'url for the TAP service')
+
+
+conf = Conf()
+
+# Now import your public class
+# Should probably have the same name as your module
+from .core import DESILegacySurvey, DESILegacySurveyClass
+
+__all__ = ['DESILegacySurvey', 'DESILegacySurveyClass',
+           'Conf', 'conf',
+           ]

astroquery/desi/core.py

@@ -0,0 +1,103 @@
+# Licensed under a 3-clause BSD style license - see LICENSE.rst
+import urllib.error
+import requests
+import warnings
+
+import pyvo as vo
+import numpy as np
+import astropy.coordinates as coord
+
+from astropy import units as u
+from astroquery.exceptions import NoResultsWarning
+from astroquery.query import BaseQuery
+from astroquery.utils import commons
+from astroquery.desi import conf
+
+__all__ = ['DESILegacySurvey', 'DESILegacySurveyClass']
+
+
+class DESILegacySurveyClass(BaseQuery):
+
+    def query_region(self, coordinates, *, width=0.5*u.arcmin, data_release=9):
+        """
+        Queries a region around the specified coordinates.
+
+        Parameters
+        ----------
+        coordinates : `~astropy.coordinates.SkyCoord`
+            coordinates around which to query.
+        width : `~astropy.units.Quantity`, optional
+            the width of the region. If missing, set to default
+            value of 0.5 arcmin.
+        data_release : int
+            the data release of the LegacySurvey to use.
+
+        Returns
+        -------
+        response : `~astropy.table.Table`
+        """
+
+        tap_service = vo.dal.TAPService(conf.tap_service_url)
+        coordinates_transformed = coordinates.transform_to(coord.ICRS)
+
+        qstr = (f"SELECT all * FROM ls_dr{data_release}.tractor WHERE "
+                f"dec<{(coordinates_transformed.dec + width).to_value(u.deg)} and "
+                f"ra>{coordinates_transformed.ra.to_value(u.deg) - width.to_value(u.deg) / np.cos(coordinates_transformed.dec)} and "
+                f"ra<{coordinates_transformed.ra.to_value(u.deg) + width.to_value(u.deg) / np.cos(coordinates_transformed.dec)}")
+
+        tap_result = tap_service.run_sync(qstr)
+        tap_result = tap_result.to_table()
+        # filter out duplicated lines from the table
+        mask = tap_result['type'] != 'D'
+        filtered_table = tap_result[mask]
+
+        return filtered_table
+
+    def get_images(self, position, *, pixels=None, width=0.5*u.arcmin, data_release=9, show_progress=True, image_band='g'):
+        """
+        Downloads the images for a certain region of interest.
+
+        Parameters
+        ----------
+        position : `astropy.coordinates`.
+            coordinates around which we define our region of interest.
+        width : `~astropy.units.Quantity`,  optional
+            the width of our region of interest.
+        data_release : int, optional
+            the data release of the LegacySurvey to use.
+        show_progress : bool, optional
+            Whether to display a progress bar if the file is downloaded
+            from a remote server.  Default is True.
+        image_band : str, optional
+            Default to 'g'
+
+        Returns
+        -------
+        list : A list of `~astropy.io.fits.HDUList` objects.
+        """
+
+        position_transformed = position.transform_to(coord.ICRS)
+
+        image_size_arcsec = width.arcsec
+        pixsize = 2 * image_size_arcsec / pixels
+
+        image_url = (f"{conf.legacysurvey_service_url}?"
+                    f"ra={position_transformed.ra.deg}&"
+                    f"dec={position_transformed.dec.deg}&"
+                    f"size={pixels}&"
+                    f"layer=ls-dr{data_release}&"
+                    f"pixscale={pixsize}&"
+                    f"bands={image_band}")
+
+        file_container = commons.FileContainer(image_url, encoding='binary', show_progress=show_progress)
+
+        try:
+            fits_file = file_container.get_fits()
+        except (requests.exceptions.HTTPError, urllib.error.HTTPError) as exp:
+            fits_file = None
+            warnings.warn(f"{str(exp)} - Problem retrieving the file at the url: {image_url}", NoResultsWarning)
+
+        return [fits_file]
+
+
+DESILegacySurvey = DESILegacySurveyClass()

astroquery/desi/tests/data/dummy_table.txt

@@ -0,0 +1,8 @@
+ra,dec,objid,type
+166.10552527002142,38.20797162140221,877,PSF
+166.10328347620825,38.211862682863625,855,PSF
+166.1146193911762,38.20826292586543,991,PSF
+166.1138080401007,38.20883307659884,3705,DUP
+166.11382707824612,38.20885008952696,982,SER
+166.11779248975387,38.211159276708706,1030,PSF
+166.11865123008005,38.2147201669633,1039,PSF

astroquery/desi/tests/setup_package.py

@@ -0,0 +1,9 @@
+# Licensed under a 3-clause BSD style license - see LICENSE.rst
+import os
+
+
+def get_package_data():
+    paths = [os.path.join('data', '*.txt'),
+             os.path.join('data', '*.fits'),
+             ]
+    return {'astroquery.desi.tests': paths}

astroquery/desi/tests/test_desi.py

@@ -0,0 +1,115 @@
+import astropy.io.votable
+import pytest
+import os
+import numpy as np
+import pyvo as vo
+
+from astropy.table import Table
+from astropy.io import fits
+from astropy import coordinates as coord
+from astroquery.utils.mocks import MockResponse
+from astroquery.utils import commons
+from astroquery import desi
+from urllib import parse
+from contextlib import contextmanager
+
+
+DATA_FILES = {
+    'dummy_tap_table': 'dummy_table.txt',
+    'dummy_tractor_fits': 'dummy_tractor.fits',
+    'dummy_hdu_list_fits': 'hdu_list_images.fits'
+}
+
+coords = coord.SkyCoord('11h04m27s +38d12m32s', frame='icrs')
+width = coord.Angle(0.5, unit='arcmin')
+pixels = 60
+data_release = 9
+emispheres_list = ['north', 'south']
+
+
+@pytest.fixture
+def patch_get(request):
+    mp = request.getfixturevalue("monkeypatch")
+
+    mp.setattr(desi.DESILegacySurvey, '_request', get_mockreturn)
+    return mp
+
+
+@pytest.fixture
+def patch_get_readable_fileobj(request):
+    @contextmanager
+    def get_readable_fileobj_mockreturn(filename, **kwargs):
+        file_obj = data_path(DATA_FILES['dummy_hdu_list_fits'])  # TODO: add images option
+        # f = None
+        with open(file_obj, 'rb') as f:
+            yield f
+
+    mp = request.getfixturevalue("monkeypatch")
+
+    mp.setattr(commons, 'get_readable_fileobj',
+               get_readable_fileobj_mockreturn)
+    return mp
+
+
+@pytest.fixture
+def patch_tap(request):
+    mp = request.getfixturevalue("monkeypatch")
+
+    mp.setattr(vo.dal.TAPService, 'run_sync', tap_mockreturn)
+    return mp
+
+
+def get_mockreturn(method, url, params=None, timeout=10, **kwargs):
+    parsed_url = parse.urlparse(url)
+    splitted_parsed_url = parsed_url.path.split('/')
+    url_filename = splitted_parsed_url[-1]
+    filename = None
+    content = None
+    if url_filename.startswith('tractor-'):
+        filename = data_path(DATA_FILES['dummy_tractor_fits'])
+
+    if filename is not None:
+        with open(filename, 'rb') as f:
+            content = f.read()
+    return MockResponse(content)
+
+
+def tap_mockreturn(url, params=None, timeout=10, **kwargs):
+    content_table = Table.read(data_path(DATA_FILES['dummy_tap_table']),
+                               format='ascii.csv', comment='#')
+    votable_table = astropy.io.votable.from_table(content_table)
+    return vo.dal.TAPResults(votable_table)
+
+
+def data_path(filename):
+    data_dir = os.path.join(os.path.dirname(__file__), 'data')
+    return os.path.join(data_dir, filename)
+
+
+def compare_result_data(result, data):
+    for col in result.colnames:
+        if result[col].dtype.type is np.string_ or result[col].dtype.type is np.str_:
+            assert np.array_equal(result[col], data[col])
+        else:
+            np.testing.assert_allclose(result[col], data[col])
+
+
+def image_tester(images, filetype):
+    assert type(images) == list
+    with fits.open(data_path(DATA_FILES[filetype])) as data:
+        assert images[0][0].header == data[0].header
+        assert np.array_equal(images[0][0].data, data[0].data)
+
+
+def test_coords_query_region(patch_tap):
+    result = desi.DESILegacySurvey.query_region(coords, width=width)
+    data = Table.read(data_path(DATA_FILES['dummy_tap_table']),
+                      format='ascii.csv', comment='#')
+    data['objid'] = data['objid'].astype(np.int64)
+    compare_result_data(result, data)
+
+
+def test_coords_get_images(patch_get_readable_fileobj, dr=data_release):
+    images_list = desi.DESILegacySurvey.get_images(coords, data_release=dr, width=width, pixels=pixels)
+
+    image_tester(images_list, 'dummy_hdu_list_fits')

astroquery/desi/tests/test_desi_remote.py

@@ -0,0 +1,41 @@
+import pytest
+
+from astroquery.desi import DESILegacySurvey
+from astropy.io.fits import HDUList
+from astropy.coordinates import SkyCoord
+from astropy.coordinates import Angle
+from astropy.table import Table
+from astroquery.exceptions import NoResultsWarning
+
+
+@pytest.mark.remote_data
+class TestLegacySurveyClass:
+
+    def test_query_region(self):
+        ra = 166.1125
+        dec = 38.209
+        coordinates = SkyCoord(ra, dec, unit='degree')
+
+        width = Angle(15, unit='arcsec')
+
+        query1 = DESILegacySurvey.query_region(coordinates, width=width, data_release=9)
+
+        assert isinstance(query1, Table)
+
+    @pytest.mark.parametrize(('ra', 'dec', 'width', 'pixels'),
+                             ((166.1125, 38.209, 0.5, 60),))
+    def test_get_images(self, ra, dec, width, pixels):
+        pos = SkyCoord(ra, dec, unit='degree')
+        width = Angle(width, unit='arcmin')
+
+        query1 = DESILegacySurvey.get_images(pos, pixels=pixels, width=width, data_release=9)
+        assert isinstance(query1, list)
+        assert isinstance(query1[0], HDUList)
+
+    def test_noresults_warning(self):
+        # Using position with no coverage
+        pos = SkyCoord(86.633212, 22.01446, unit='degree')
+        width = Angle(3, unit='arcmin')
+
+        with pytest.warns(NoResultsWarning):
+            DESILegacySurvey.get_images(pos, width=width, pixels=100)

astroquery/utils/commons.py

@@ -24,10 +24,8 @@
 from ..exceptions import TimeoutError, InputWarning
 from .. import version
 
-
 CoordClasses = (SkyCoord, BaseCoordinateFrame)
 
-
 __all__ = ['send_request',
            'parse_coordinates',
            'TableList',