diff --git a/tools/mwa/.shed.yml b/tools/mwa/.shed.yml new file mode 100644 index 00000000..0ff4f255 --- /dev/null +++ b/tools/mwa/.shed.yml @@ -0,0 +1,9 @@ +categories: +- Astronomy +description: MWA +homepage_url: null +long_description: MWA +name: mwa_astro_tool +owner: astroteam +remote_repository_url: https://github.com/esg-epfl-apc/tools-astro/tree/main/tools +type: unrestricted diff --git a/tools/mwa/Image.py b/tools/mwa/Image.py new file mode 100644 index 00000000..d27506da --- /dev/null +++ b/tools/mwa/Image.py @@ -0,0 +1,91 @@ +#!/usr/bin/env python +# coding: utf-8 + +# flake8: noqa + +import json +import os +import shutil + +import astropy.units as u +import matplotlib.pyplot as plt +from astropy.wcs import WCS +from astroquery.skyview import SkyView +from oda_api.data_products import ImageDataProduct, PictureProduct +from oda_api.json import CustomJSONEncoder + +src_name = "1ES 0229+200" # http://odahub.io/ontology#AstrophysicalObject +RA = 38.202562 # http://odahub.io/ontology#PointOfInterestRA +DEC = 20.288191 # http://odahub.io/ontology#PointOfInterestDEC +T1 = "2000-10-09T13:16:00.0" # http://odahub.io/ontology#StartTime +T2 = "2022-10-10T13:16:00.0" # http://odahub.io/ontology#EndTime +Radius = 1.0 # http://odahub.io/ontology#AngleDegrees +pixsize = 0.01 # http://odahub.io/ontology#AngleDegrees +Frequency = "GLEAM 170-231 MHz" # http://odahub.io/ontology#String ; oda:allowed_value "GLEAM 72-103 MHz","GLEAM 103-134 MHz","GLEAM 139-170 MHz","GLEAM 170-231 MHz" + +_galaxy_wd = os.getcwd() + +with open("inputs.json", "r") as fd: + inp_dic = json.load(fd) +if "_data_product" in inp_dic.keys(): + inp_pdic = inp_dic["_data_product"] +else: + inp_pdic = inp_dic + +for vn, vv in inp_pdic.items(): + if vn != "_selector": + globals()[vn] = type(globals()[vn])(vv) + +pixels = int(2 * Radius / pixsize) + 1 +Radius *= u.deg +pos = str(RA) + ", " + str(DEC) +pixels + +hdul = SkyView.get_images( + position=pos, survey=[Frequency], pixels=pixels, radius=Radius +) + +hdu = hdul[0] +hdu[0].header +wcs = WCS(hdu[0].header) + +image = hdu[0].data + +ax = plt.subplot(projection=wcs) +im = ax.imshow(image, origin="lower") +ax.coords.grid(True, color="white", ls="solid") +plt.colorbar(im, label="Jy/beam") +plt.savefig("Image.png", format="png", bbox_inches="tight") + +hdu.writeto("Image.fits", overwrite=True) +bin_image = PictureProduct.from_file("Image.png") +fits_image = ImageDataProduct.from_fits_file("Image.fits") + +picture = bin_image # http://odahub.io/ontology#ODAPictureProduct +image = fits_image # http://odahub.io/ontology#Image + +# output gathering +_galaxy_meta_data = {} +_oda_outs = [] +_oda_outs.append(("out_Image_picture", "picture_galaxy.output", picture)) +_oda_outs.append(("out_Image_image", "image_galaxy.output", image)) + +for _outn, _outfn, _outv in _oda_outs: + _galaxy_outfile_name = os.path.join(_galaxy_wd, _outfn) + if isinstance(_outv, str) and os.path.isfile(_outv): + shutil.move(_outv, _galaxy_outfile_name) + _galaxy_meta_data[_outn] = {"ext": "_sniff_"} + elif getattr(_outv, "write_fits_file", None): + _outv.write_fits_file(_galaxy_outfile_name) + _galaxy_meta_data[_outn] = {"ext": "fits"} + elif getattr(_outv, "write_file", None): + _outv.write_file(_galaxy_outfile_name) + _galaxy_meta_data[_outn] = {"ext": "_sniff_"} + else: + with open(_galaxy_outfile_name, "w") as fd: + json.dump(_outv, fd, cls=CustomJSONEncoder) + _galaxy_meta_data[_outn] = {"ext": "json"} + +with open(os.path.join(_galaxy_wd, "galaxy.json"), "w") as fd: + json.dump(_galaxy_meta_data, fd) +print("*** Job finished successfully ***") diff --git a/tools/mwa/Spectrum.py b/tools/mwa/Spectrum.py new file mode 100644 index 00000000..c1c59b20 --- /dev/null +++ b/tools/mwa/Spectrum.py @@ -0,0 +1,117 @@ +#!/usr/bin/env python +# coding: utf-8 + +# flake8: noqa + +import json +import os +import shutil + +import astropy.units as u +import matplotlib.pyplot as plt +import numpy as np +from astropy.constants import h +from oda_api.data_products import ODAAstropyTable, PictureProduct +from oda_api.json import CustomJSONEncoder +from pyvo import registry # version >=1.4.1 + +src_name = "1ES 0229+200" # http://odahub.io/ontology#AstrophysicalObject +RA = 38.202562 # http://odahub.io/ontology#PointOfInterestRA +DEC = 20.288191 # http://odahub.io/ontology#PointOfInterestDEC +T1 = "2000-10-09T13:16:00.0" # http://odahub.io/ontology#StartTime +T2 = "2022-10-10T13:16:00.0" # http://odahub.io/ontology#EndTime +Radius = 0.05 # http://odahub.io/ontology#AngleDegrees + +_galaxy_wd = os.getcwd() + +with open("inputs.json", "r") as fd: + inp_dic = json.load(fd) +if "_data_product" in inp_dic.keys(): + inp_pdic = inp_dic["_data_product"] +else: + inp_pdic = inp_dic + +for vn, vv in inp_pdic.items(): + if vn != "_selector": + globals()[vn] = type(globals()[vn])(vv) + +conesearch_radius = Radius # in degrees +conesearch_center = (RA, DEC) + +# the catalogue name in VizieR +CATALOGUE = "VIII/100" +# each resource in the VO has an identifier, called ivoid. For vizier catalogs, +# the VO ids can be constructed like this: +catalogue_ivoid = f"ivo://CDS.VizieR/{CATALOGUE}" +# the actual query to the registry +voresource = registry.search(ivoid=catalogue_ivoid)[0] + +conesearch_records = voresource.get_service("conesearch").search( + pos=conesearch_center, + sr=conesearch_radius, +) +conesearch_records + +h_p = (h / u.s).to(u.eV).value # Planck constant in eV*s + +conesearch_records.fieldnames +nu = [] +F = [] +F_err = [] +for f in conesearch_records.fieldnames: + if (f[:4] == "Fint") and (f[4] != "w") and (f[4] != "f"): + nu.append(int(f[-3:]) * 1e6) # in Hz + F.append(conesearch_records[f][0] * nu[-1] * 1e-23) + F_err.append(conesearch_records["e_" + f][0] * nu[-1] * 1e-23) + +E = h_p * np.array(nu) +plt.errorbar(E, F, F_err) +plt.xscale("log") +plt.yscale("log") +plt.xlabel("$E$, eV") +plt.ylabel("$E F_E$, erg/cm$^2$s") +plt.savefig("Spectrum.png", format="png", bbox_inches="tight") + +bin_image = PictureProduct.from_file("Spectrum.png") +from astropy.table import Table + +data = [E, F, F_err] +names = ("E[eV]", "Flux[erg/cm2s]", "Flux_error[erg/cm2s]") +spec = ODAAstropyTable(Table(data, names=names)) + +picture_png = bin_image # http://odahub.io/ontology#ODAPictureProduct +spectrum_astropy_table = spec # http://odahub.io/ontology#ODAAstropyTable + +# output gathering +_galaxy_meta_data = {} +_oda_outs = [] +_oda_outs.append( + ("out_Spectrum_picture_png", "picture_png_galaxy.output", picture_png) +) +_oda_outs.append( + ( + "out_Spectrum_spectrum_astropy_table", + "spectrum_astropy_table_galaxy.output", + spectrum_astropy_table, + ) +) + +for _outn, _outfn, _outv in _oda_outs: + _galaxy_outfile_name = os.path.join(_galaxy_wd, _outfn) + if isinstance(_outv, str) and os.path.isfile(_outv): + shutil.move(_outv, _galaxy_outfile_name) + _galaxy_meta_data[_outn] = {"ext": "_sniff_"} + elif getattr(_outv, "write_fits_file", None): + _outv.write_fits_file(_galaxy_outfile_name) + _galaxy_meta_data[_outn] = {"ext": "fits"} + elif getattr(_outv, "write_file", None): + _outv.write_file(_galaxy_outfile_name) + _galaxy_meta_data[_outn] = {"ext": "_sniff_"} + else: + with open(_galaxy_outfile_name, "w") as fd: + json.dump(_outv, fd, cls=CustomJSONEncoder) + _galaxy_meta_data[_outn] = {"ext": "json"} + +with open(os.path.join(_galaxy_wd, "galaxy.json"), "w") as fd: + json.dump(_galaxy_meta_data, fd) +print("*** Job finished successfully ***") diff --git a/tools/mwa/mwa_astro_tool.xml b/tools/mwa/mwa_astro_tool.xml new file mode 100644 index 00000000..cf32759e --- /dev/null +++ b/tools/mwa/mwa_astro_tool.xml @@ -0,0 +1,102 @@ + + + ipython + astropy + gammapy + matplotlib + astroquery + scipy + oda-api + nbconvert + wget + pandas + + + python '$__tool_directory__/${_data_product._selector}.py' + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + _data_product['_selector'] == 'Image' + + + _data_product['_selector'] == 'Image' + + + _data_product['_selector'] == 'Spectrum' + + + _data_product['_selector'] == 'Spectrum' + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + This is help for MWA + + + 10.1093/mnras/stw2337 + + \ No newline at end of file