{ "cells": [ { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "# Astropy and the VO" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "This is supposed to be an appetiser (and, to some extent, quarry for sample code) for doing\n", "fun and profitable stuff with VAO’s astropy extension for VO access (called pyVO) as well as\n", "astropy’s built-in VO capabilities.\n", "\n", "PyVO is available at https://pyvo.readthedocs.io/en/latest/index.html#installation." ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "## 1. Query a single SIAP service" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "SIAP is the VO’s protocol to access image servers. Here’s how you’d do a search at multiple\n", "positions with date constraints:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "from astropy.time import Time\n", "from pyvo.dal import sia\n", "import warnings\n", "warnings.filterwarnings('ignore', module=\"astropy.io.votable.*\")\n", "warnings.filterwarnings('ignore', module=\"astropy.units.format.utils\")" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "ACCESS_URL = \"http://dc.g-vo.org/lswscans/res/positions/siap/siap.xml?\"\n", "DATE_MIN = Time(\"1902-01-01\", scale=\"tt\").mjd\n", "DATE_MAX = Time(\"1922-12-31T23:59:59\", scale=\"tt\").mjd\n", "\n", "svc = sia.SIAService(ACCESS_URL)\n", "for pos in [\n", " (10, 20),\n", " (240, -10),\n", " (45, 85)]:\n", " images = svc.search(pos, (0.5, 0.5), verbosity=2) #returns VOTable \n", " \n", " print(\"Available images and respective dates:\")\n", " for rec_no, match in enumerate(images):\n", " \n", " if not DATE_MINmax_queries:\n", " break\n", " try:\n", " search_one_service(svc.access_url)\n", " except:\n", " import traceback; traceback.print_exc()" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "Note: The exception catcher is there since not all services claiming to be standards-compliant actually\n", "are. It doesn’t hurt to complain to the service operators if a service you’re interested in behaves\n", "weirdly – sometimes the operators haven’t noticed it’s broken or it just broke." ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "## 3. Add SAMP Magic" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "Let’s send over images we’re interested in via SAMP.\n", "\n", "SAMP-enabling programs may not come quite natural to people that so far have mainly written\n", "fairly linear science code. But here it’s no big deal except for the necessity to manage the\n", "connection to the SAMP hub (with e.g. TOPCAT), which is taken care of by a hand-crafted context manager in our\n", "code.\n", "Also, the way arguments are passed between SAMP services may seem a bit funky – but think\n", "of the mtype as the function name, and passing arguments in dictionaries instead of tuples isn’t\n", "that far-fetched, either." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "from contextlib import contextmanager\n", "from astropy.samp import SAMPIntegratedClient" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "\"\"\"\n", "A little script doing an all-VO SIAP query for some positions and a date \n", "range; the results can be sent to SAMP clients.\n", "\"\"\"\n", "\n", "max_queries = 20 # desired number of queries\n", "\n", "DATE_MIN = Time(\"1990-01-01\", scale=\"tt\").mjd\n", "DATE_MAX = Time(\"2005-12-31T23:59:59\", scale=\"tt\").mjd\n", "NAME_TEMPLATE = \"global-{:03}\"\n", "\n", "class ImageSender(object):\n", " \"\"\"A SAMP handler broadcasting images.\n", "\n", " Construct this class with a SAMP connection and call broadcast_image\n", " to have SAMP clients load the images.\n", " \"\"\"\n", " def __init__(self, client):\n", " self.client = client\n", " self.image_count = 0\n", "\n", " def broadcast_image(self, image_url):\n", " message = {\n", " \"samp.mtype\": \"table.load.fits\",\n", " \"samp.params\": {\n", " \"url\": image_url,\n", " \"name\": NAME_TEMPLATE.format(self.image_count),\n", " }}\n", " self.client.notify_all(message)\n", "\n", "\n", "def search_one_service(access_url, image_sender):\n", " print(\"Now querying \", access_url)\n", " svc = sia.SIAService(access_url)\n", " for pos in [\n", " (10, 20)]:\n", " images = svc.search(pos, (0.5, 0.5), verbosity=2)\n", " dateName = images.fieldname_with_ucd(\"VOX:Image_MJDateObs\")\n", " if dateName is None:\n", " return\n", "\n", " for rec_no, match in enumerate(images):\n", " \n", " if not DATE_MINmax_queries:\n", " break\n", " try:\n", " search_one_service(svc.access_url, image_sender)\n", " except Exception:\n", " import traceback; traceback.print_exc()\n", "\n" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "## 4. Higher Magic" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "Let’s say you’re debugging your pipeline and want to manually inspect “weird” objects by\n", "checking what a set of other catalogs have on them.\n", "\n", "So, a script would have to notice when a table row is selected:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "class VicinitySearcher(object):\n", " def __init__(self, client):\n", " self.client.bind_receive_call(\n", " \"table.load.votable\", self.load_VOTable)\n", " self.client.bind_receive_notification(\"table.highlight.row\",\n", " self.handle_selection)\n", " def load_VOTable(self, private_key, ...\n", " self.cur_table = Table.read(params[’url’])\n", " self.cur_id = params[\"table-id\"]\n", " self.client.reply(msg_id,\n", " {\"samp.status\": \"samp.ok\", \"samp.result\": {}})\n", " def handle_selection(self, private_key, ...\n", " if params[\"table-id\"]!=self.cur_id:\n", " return\n", " table_index = int(params[\"row\"])\n", " response = self.make_response_table(table_index)\n", " if response is not None:\n", " self.send_table(response)" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "Essentially, we’re telling the SAMP hub we’re interested in highlights of rows in tables. However,\n", "as multiple tables may be around, in our handler we make sure we only do something if we’re\n", "dealing with the last table we got sent." ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "## 5. United Services of the VO" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "If we query multiple services and want to get a single table back, we must combine their results.\n", "The Cone Search standard (see https://www.ivoa.net/documents/REC/DAL/ConeSearch-20080222.html) and UCDs (see http://docs.g-vo.org/pyvo/html/page014.html) help:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "def make_response_table(self, table_index):\n", " ras, decs, pmras, pmdecs, svcs = [], [], [], [], []\n", " for service in self.services:\n", " cone_result = service.search((ra,dec),\n", " self.vicinity_size/3600.)\n", " ras.extend(cone_result.getcolumn(\n", " cone_result.fieldname_with_ucd(\"POS_EQ_RA_MAIN\")))\n", " decs.extend(cone_result.getcolumn(\n", " cone_result.fieldname_with_ucd(\"POS_EQ_DEC_MAIN\")))\n", " pmra_name = cone_result.fieldname_with_ucd(\"pos.pm;pos.eq.ra\")\n", " if pmra_name:\n", " pmras.extend(cone_result.getcolumn(pmra_name))\n", " else:\n", " pmras.extend([None]*nrecs)" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "The UCDs (those for RA and DEC look odd because the cone search standard is so ancient)\n", "help us retrieve columns by physics (what’s the RA? what’s the PM in RA?).\n", "Incidentally, the units on the PMs are not defined by the standard, and to actually unify them\n", "we’d need to do unit conversion. That’s standard astropy, though.\n", "So, finally, here’s the full program with a usage scenario in the module docstring:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import os\n", "import tempfile\n", "import traceback\n", "\n", "from astropy.table import Table\n", "from pyvo.dal import scs" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "\"\"\"\n", "A quick example showing astropy and pyvo working hand in hand with the\n", "rest of the VO.\n", "\n", "This program expects Aladin to run. It then waits for a tables to be sent,\n", "and when a row is selected, it will search some (SERVICE_META) cone\n", "search services. The results are joined and sent to aladin with\n", "positions, proper motions, and source.\n", "\n", "Sample use:\n", "\n", "(1) start TOPCAT, Aladin and run this cell\n", "(2) in TOPCAT, open VO/Cone Search, look for \"transitional YSOs\" as Keyword\n", "(3) select the Magnier+ 1999 service (Sort Name: J/A+A/352/228), go the Cone \n", " Parameters and make RA and Dec 0 degrees, Radius 180 degrees and hit \"ok\"\n", "(4) broadcast table in the TOPCAT's main window (3rd button from the left)\n", "(5) in Aladin, pan and zoom until you have a catalog object in a FoV of\n", " an arcminute or so\n", "(6) hover over the object to pull in the potential matches\n", "(7) select the items to see the catalog entries.\n", "\"\"\"\n", "\n", "\n", "\n", "SERVICE_META = [\n", " (\"PPMXL\", \"http://dc.zah.uni-heidelberg.de/ppmxl/q/cone/scs.xml?\"),\n", " (\"2MASS\", \"http://dc.zah.uni-heidelberg.de/2mass/res/2mass/q/scs.xml?\"),\n", " (\"UCAC4\", \"http://dc.zah.uni-heidelberg.de/ucac4/q/s/scs.xml?\")]\n", "\n", "\n", "\n", "@contextmanager\n", "def samp_accessible(astropy_table):\n", " \"\"\"a context manager making astroy_table available under a (file)\n", " URL for the controlled section.\n", "\n", " This is useful with uploads.\n", " \"\"\"\n", " handle, f_name = tempfile.mkstemp(suffix=\".xml\")\n", " with os.fdopen(handle, \"w\") as f:\n", " astropy_table.write(output=f,\n", " format=\"votable\")\n", " try:\n", " yield \"file://\"+f_name\n", " finally:\n", " os.unlink(f_name)\n", "\n", "\n", "def get_name_for_ucd(ucd, table):\n", " \"\"\"returns the name of a column having ucd in table.\n", "\n", " If there are multiple such columns, a random one is returned. If there\n", " are none, a key error is raised.\n", " \"\"\"\n", " ucd = ucd.lower()\n", "\n", " for col in table.columns.values():\n", " if col.meta.get(\"ucd\").lower()==ucd:\n", " return col.name\n", " raise KeyError(ucd)\n", "\n", "\n", "class VicinitySearcher(object):\n", " \"\"\"The SAMP handling class.\n", "\n", " This is where the action takes place: receiving VOTables, handling\n", " notifications of selected rows, querying the remote services.\n", "\n", " True, in a less one-off program this should be less god-like, and\n", " at least make_response_table shouln't be part of this.\n", " \"\"\"\n", " vicinity_size = 30 # arcsec\n", " client_name = \"Aladin\" # samp.name of the client for the match table\n", "\n", " def __init__(self, client):\n", " self.client = client\n", " self.cur_table = self.cur_id = None\n", " self.dest_client = self.find_client(self.client_name)\n", "\n", " self.services = []\n", " for short_name, access_url in SERVICE_META:\n", " self.services.append(scs.SCSService(access_url))\n", " self.services[-1].my_tag = short_name\n", "\n", " self.client.bind_receive_call(\n", " \"table.load.votable\", self.load_VOTable)\n", " self.client.bind_receive_notification(\"table.highlight.row\", \n", " self.handle_selection)\n", "\n", " def load_VOTable(self, private_key, sender_id, msg_id, mtype, \n", " params, extra):\n", " \"\"\"the SAMP handler to load VOTables.\n", "\n", " (binding is done in the constructor)\n", " \"\"\"\n", " self.cur_table = Table.read(params['url'])\n", " self.ra_name = get_name_for_ucd(\"POS_EQ_RA_MAIN\", self.cur_table)\n", " self.dec_name = get_name_for_ucd(\"POS_EQ_DEC_MAIN\", self.cur_table)\n", " self.cur_id = params[\"table-id\"]\n", "\n", " self.client.reply(msg_id, \n", " {\"samp.status\": \"samp.ok\", \"samp.result\": {}})\n", "\n", " def handle_selection(self, private_key, sender_id, msg_id, mtype, \n", " params, extra):\n", " \"\"\"the SAMP handler for a row selection in our current table.\n", " \"\"\"\n", " with debug():\n", " if params[\"table-id\"]!=self.cur_id:\n", " # doesn't concern us, not our table\n", " return\n", " table_index = int(params[\"row\"])\n", " print(\"Row selected:\", table_index)\n", " response = self.make_response_table(table_index)\n", " \n", " if response is not None:\n", " self.send_table(response)\n", " \n", " def send_table(self, response_table):\n", " \"\"\"sends the (astropy) response_table via SAMP.\n", " \"\"\"\n", " with samp_accessible(response_table) as table_url:\n", " message = {\n", " \"samp.mtype\": \"table.load.votable\",\n", " \"samp.params\": {\n", " \"url\": table_url,\n", " }}\n", " self.client.call_and_wait(self.dest_client, message, \"10\")\n", "\n", " def find_client(self, samp_name):\n", " \"\"\"returns the SAMP id of the client with samp.name samp_name.\n", "\n", " This will raise a KeyError if the client is not on the hub.\n", " \"\"\"\n", " for client_id in self.client.get_registered_clients():\n", " if self.client.get_metadata(client_id).get(\"samp.name\")==samp_name:\n", " return client_id\n", " raise KeyError(samp_name)\n", " \n", " def make_response_table(self, table_index):\n", " \"\"\"returns an astropy table (or None) for the row table_index.\n", "\n", " This is essentially the \"user code\" that reacts on the incoming\n", " messages.\n", " \"\"\"\n", " ra = self.cur_table[self.ra_name][table_index]\n", " dec = self.cur_table[self.dec_name][table_index]\n", "\n", " ras, decs, pmras, pmdecs, svcs = [], [], [], [], []\n", " for service in self.services:\n", " print (\"Querying \", service.my_tag)\n", " cone_result = service.search((ra,dec),\n", " self.vicinity_size/3600.)\n", " # len(core_result) appears broken ATM\n", " nrecs = len(cone_result.getcolumn(\n", " cone_result.fieldname_with_ucd(\"POS_EQ_RA_MAIN\")))\n", "\n", " ras.extend(cone_result.getcolumn(\n", " cone_result.fieldname_with_ucd(\"POS_EQ_RA_MAIN\")))\n", " decs.extend(cone_result.getcolumn(\n", " cone_result.fieldname_with_ucd(\"POS_EQ_DEC_MAIN\")))\n", "\n", " pmra_name = cone_result.fieldname_with_ucd(\"pos.pm;pos.eq.ra\")\n", " if pmra_name:\n", " # TODO: check unit and covert if necessary\n", " pmras.extend(cone_result.getcolumn(pmra_name))\n", " else:\n", " pmras.extend([None]*nrecs)\n", "\n", " pmdec_name = cone_result.fieldname_with_ucd(\"pos.pm;pos.eq.dec\")\n", " if pmdec_name:\n", " pmdecs.extend(cone_result.getcolumn(pmdec_name))\n", " else:\n", " pmdecs.extend([None]*nrecs)\n", " \n", " svcs.extend([service.my_tag]*nrecs)\n", "\n", " if not ras:\n", " return None\n", " else:\n", " print (\"Found {} matches\".format(len(ras)))\n", "\n", " res = Table([ras, decs, pmras, pmdecs, svcs],\n", " names=(\"ra\", \"dec\", \"pmra\", \"pmdec\", \"svc\"))\n", " res[\"ra\"].unit = \"deg\"\n", " res[\"dec\"].unit = \"deg\"\n", " res[\"pmra\"].unit = \"deg/yr\"\n", " res[\"pmdec\"].unit = \"deg/yr\"\n", " # TODO: enter UCDs in the meta dict\n", "\n", " return res\n", "\n", "\n", "@contextmanager\n", "def SAMP_conn():\n", " \"\"\"a context manager to give the controlled block a SAMP connection.\n", "\n", " The program will disconnect as the controlled block is exited.\n", " \"\"\"\n", " client = SAMPIntegratedClient(name=\"smpsmp\", \n", " description=\"A sample for the use of SAMP\")\n", " client.connect()\n", " try:\n", " yield client\n", " finally:\n", " client.disconnect()\n", "\n", "\n", "def main():\n", " with SAMP_conn() as client:\n", " handler = VicinitySearcher(client)\n", " print (\"Listening. Send me a table, hit return to exit.\")\n", " input()\n", " \n", "\n", "if __name__==\"__main__\":\n", " main()" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "The application of the little program above should look along the lines of the following screenshot:" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "The following method can optionally be included if tracebacks of exceptions should be displayed." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "@contextmanager\n", "def debug():\n", " \"\"\"dumps tracebacks of exceptions within the controlled block.\n", "\n", " That's useful here as astropy's SAMP silently swallows exceptions\n", " in handlers.\n", " \"\"\"\n", " try:\n", " yield\n", " except:\n", " traceback.print_exc()\n", " raise" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "When riding home. . ." ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ ". . . take this with you: http://docs.g-vo.org/pyvo.pdf." ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3" }, "vscode": { "interpreter": { "hash": "916dbcbb3f70747c44a77c7bcd40155683ae19c65e1c03b4aa3499c5328201f1" } } }, "nbformat": 4, "nbformat_minor": 2 }