Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Update tool HESS to 0.0.2+galaxy0 #94

Merged
merged 24 commits into from
Apr 18, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension


Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
25 changes: 24 additions & 1 deletion .github/workflows/live-preview.yaml
Original file line number Diff line number Diff line change
Expand Up @@ -7,6 +7,8 @@ on:
- opened
- labeled
- synchronize
- unlabeled
- closed
paths-ignore:
- '.github/**'
- 'deprecated/**'
Expand All @@ -18,7 +20,7 @@ env:

jobs:
setup:
if: contains( github.event.pull_request.labels.*.name, 'test-live')
if: contains( github.event.pull_request.labels.*.name, 'test-live') && github.event.action != 'closed'
name: Determine changed repositories
runs-on: ubuntu-latest
outputs:
Expand Down Expand Up @@ -73,3 +75,24 @@ jobs:

ansible-playbook deploy.yml
done

uninstall:
name: Uninstall
if: contains( github.event.pull_request.labels.*.name, 'test-live') == false || github.event.action == 'closed'
runs-on: ubuntu-latest

steps:
- uses: actions/checkout@v3
- name: Install requirenments
run: sudo apt-get install -y ansible
- name: Install SSH key
uses: shimataro/ssh-key-action@v2
with:
key: ${{ secrets.SSH_PRIVATE_KEY }}
known_hosts: ${{ secrets.KNOWN_HOSTS }}
- name: Uninstall tools
env:
pr_num: ${{ github.event.pull_request.number }}
run: |
cd deploy-preview
ansible-playbook remove.yml
1 change: 0 additions & 1 deletion deploy-preview/deploy.yml
Original file line number Diff line number Diff line change
Expand Up @@ -44,7 +44,6 @@
path: "{{ dest_tools_conf }}"
regexp: '^\s*<tool\s+file="{{ tool_dir }}_pr{{ pr_num }}/'
state: absent
loop: "{{ tools_list }}"

- name: Enable tools
lineinfile:
Expand Down
26 changes: 26 additions & 0 deletions deploy-preview/remove.yml
Original file line number Diff line number Diff line change
@@ -0,0 +1,26 @@
---
- hosts: all
become: true
become_user: root
vars:
dest_tools_dir: /srv/galaxy/staging_tools
dest_tools_conf: /srv/galaxy/server/config/staging_tool_conf.xml
pr_num: "{{ lookup('env', 'pr_num') }}"


tasks:
- name: Clean reference
lineinfile:
path: "{{ dest_tools_conf }}"
regexp: '^\s*<tool\s+file=".*_pr{{ pr_num }}/'
state: absent

- name: Restart Galaxy
command: galaxyctl restart
changed_when: true

- name: Remove tool dir
ansible.builtin.file:
path: "{{ dest_tools_dir }}/.*_pr{{ pr_num }}"
state: absent
mode: 0755
83 changes: 53 additions & 30 deletions tools/hess/Image.py
Original file line number Diff line number Diff line change
Expand Up @@ -14,9 +14,13 @@
from astropy.io import fits
from astropy.time import Time
from numpy import cos, pi
from oda_api.api import ProgressReporter
from oda_api.data_products import ImageDataProduct, PictureProduct
from oda_api.json import CustomJSONEncoder

pr = ProgressReporter()
pr.report_progress(stage="Progress", progress=5.0)

if os.path.exists("hess_dl3_dr1.tar.gz") == False:
get_ipython().system( # noqa: F821
"wget https://zenodo.org/record/1421099/files/hess_dl3_dr1.tar.gz"
Expand All @@ -28,12 +32,12 @@
DEC = 22.014700 # 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 = 2.5 # http://odahub.io/ontology#AngleDegrees
Radius = 1.0 # http://odahub.io/ontology#AngleDegrees
pixsize = (
0.1 # http://odahub.io/ontology#AngleDegrees ; oda:label "Pixel size"
0.05 # http://odahub.io/ontology#AngleDegrees ; oda:label "Pixel size"
)
Emin = 100.0 # http://odahub.io/ontology#Energy_GeV
Emax = 10000.0 # http://odahub.io/ontology#Energy_GeV
Emin = 1 # http://odahub.io/ontology#Energy_TeV
Emax = 100.0 # http://odahub.io/ontology#Energy_TeV

_galaxy_wd = os.getcwd()

Expand Down Expand Up @@ -95,9 +99,14 @@
message

cdec = cos(DEC * pi / 180.0)
Npix = int(4 * Radius / pixsize) + 1
RA_bins = np.linspace(RA - Radius / cdec, RA + Radius / cdec, Npix + 1)
DEC_bins = np.linspace(DEC - Radius, DEC + Radius, Npix + 1)
Npix = int(2 * Radius / pixsize) + 1
RA_bins = np.linspace(
RA - Npix * pixsize / cdec / 2, RA + Npix * pixsize / cdec / 2, Npix + 1
)
DEC_bins = np.linspace(
DEC - Npix * pixsize / 2, DEC + Npix * pixsize / 2, Npix + 1
)

image = np.zeros((Npix, Npix))
for f in OBSlist:
hdul = fits.open("data/" + f)
Expand All @@ -106,62 +115,76 @@
ev_dec = ev["DEC"]
ev_en = ev["ENERGY"]
ev_time = ev["TIME"]
h = np.histogram2d(ev_ra, ev_dec, bins=[RA_bins, DEC_bins])
mask = (ev_en > Emin) & (ev_en < Emax)
h = np.histogram2d(ev_ra[mask], ev_dec[mask], bins=[RA_bins, DEC_bins])
image += h[0]
hdul.close()

image = np.transpose(image)

plt.imshow(
np.flip(image, axis=1),
extent=(RA_bins[-1], RA_bins[0], DEC_bins[0], DEC_bins[-1]),
image,
extent=(RA_bins[0], RA_bins[-1], DEC_bins[0], DEC_bins[-1]),
origin="lower",
)
plt.colorbar()

plt.xlim(*plt.xlim()[::-1])

plt.xlabel("RA, degrees")
plt.ylabel("DEC,degrees")
plt.savefig("Image.png", format="png")

# Create a new WCS object. The number of axes must be set
# from the start
w = wcs.WCS(naxis=2)

# Set up an "Airy's zenithal" projection
# Vector properties may be set with Python lists, or Numpy arrays
w.wcs.crpix = [Npix / 2.0, Npix / 2.0]
w.wcs.cdelt = np.array([pixsize / cdec, pixsize])
w.wcs.crval = [RA, DEC]
w.wcs.ctype = ["RA---AIR", "DEC--AIR"]
w.wcs.set_pv([(2, 1, 45.0)])
w.wcs.ctype = ["RA---CAR", "DEC--CAR"]
# we need a Plate carrée (CAR) projection since histogram is binned by ra-dec
# the peculiarity here is that CAR projection produces rectilinear grid only if CRVAL2==0
# also, we will follow convention of RA increasing from right to left (CDELT1<0, need to flip an input image)
# otherwise, aladin-lite doesn't show it
w.wcs.crval = [RA, 0]
w.wcs.crpix = [Npix / 2.0 + 0.5, 0.5 - DEC_bins[0] / pixsize]
w.wcs.cdelt = np.array([-pixsize / cdec, pixsize])

# Now, write out the WCS object as a FITS header
header = w.to_header()

# header is an astropy.io.fits.Header object. We can use it to create a new
# PrimaryHDU and write it to a file.
hdu = fits.PrimaryHDU(image, header=header)
hdu = fits.PrimaryHDU(np.flip(image, axis=1), header=header)
hdu.writeto("Image.fits", overwrite=True)
hdu = fits.open("Image.fits")
im = hdu[0].data
from astropy.wcs import WCS

wcs = WCS(hdu[0].header)
plt.subplot(projection=wcs)
wcs1 = wcs.WCS(hdu[0].header)
ax = plt.subplot(projection=wcs1)
plt.imshow(im, origin="lower")
plt.colorbar(label="Counts per pixel")
plt.scatter(
[RA], [DEC], marker="x", color="white", transform=ax.get_transform("world")
)
plt.text(
RA,
DEC + 0.5 * pixsize,
src_name,
color="white",
transform=ax.get_transform("world"),
)

plt.grid(color="white", ls="solid")
plt.xlabel("RA")
plt.ylabel("Dec")
pr.report_progress(stage="Progress", progress=100.0)
plt.savefig("Image.png", format="png")

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
png = bin_image # http://odahub.io/ontology#ODAPictureProduct
fits = 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))
_oda_outs.append(("out_Image_png", "png_galaxy.output", png))
_oda_outs.append(("out_Image_fits", "fits_galaxy.output", fits))

for _outn, _outfn, _outv in _oda_outs:
_galaxy_outfile_name = os.path.join(_galaxy_wd, _outfn)
Expand Down
69 changes: 43 additions & 26 deletions tools/hess/Lightcurve.py
Original file line number Diff line number Diff line change
Expand Up @@ -21,17 +21,21 @@
"wget https://zenodo.org/record/1421099/files/hess_dl3_dr1.tar.gz"
)
get_ipython().system("tar -zxvf hess_dl3_dr1.tar.gz") # noqa: F821
from oda_api.api import ProgressReporter

pr = ProgressReporter()
pr.report_progress(stage="Progress", progress=5.0)

src_name = "Crab" # http://odahub.io/ontology#AstrophysicalObject
RA = 83.628700 # http://odahub.io/ontology#PointOfInterestRA
DEC = 22.014700 # http://odahub.io/ontology#PointOfInterestDEC
T1 = "2003-10-09T13:16:00.0" # http://odahub.io/ontology#StartTime
T2 = "2005-10-10T13:16:00.0" # http://odahub.io/ontology#EndTime
T1 = "2004-11-20T13:16:00.0" # http://odahub.io/ontology#StartTime
T2 = "2004-12-20T13:16:00.0" # http://odahub.io/ontology#EndTime
Radius = 2.5 # http://odahub.io/ontology#AngleDegrees
R_s = 0.2 # http://odahub.io/ontology#AngleDegrees
Emin = 100.0 # http://odahub.io/ontology#Energy_GeV
Emax = 10000.0 # http://odahub.io/ontology#Energy_GeV
NTbins = 10 # http://odahub.io/ontology#Integer
Emin = 1 # http://odahub.io/ontology#Energy_TeV
Emax = 100.0 # http://odahub.io/ontology#Energy_TeV
NTbins = 30 # http://odahub.io/ontology#Integer

_galaxy_wd = os.getcwd()

Expand Down Expand Up @@ -105,7 +109,12 @@
flux_b = np.zeros(NTbins)
flux_b_err = np.zeros(NTbins)
Expos = np.zeros(NTbins)
src_cts = np.zeros(NTbins)
bkg_cts = np.zeros(NTbins)
for count, f in enumerate(OBSlist):
pr.report_progress(
stage="Progress", progress=5.0 + 95 * count / len(OBSlist)
)
hdul = fits.open("data/" + f)
RA_pnt = hdul[1].header["RA_PNT"]
DEC_pnt = hdul[1].header["DEC_PNT"]
Expand All @@ -124,6 +133,7 @@
ev_dec = ev["DEC"]
ev_en = ev["ENERGY"]
ev_time = (ev["TIME"] - Trun_start) / 86400.0 + Tbegs[count]
Nevents = len(ev)
print(ev_time[0])
ev_coords = SkyCoord(ev_ra, ev_dec, unit="degree")
sep_s = ev_coords.separation(Coords_s).deg
Expand All @@ -138,18 +148,31 @@
Thmin = hdu["THETA_LO"][0]
Thmax = hdu["THETA_HI"][0]
ind = np.argmin((Thmin - dist) ** 2)
AA = AA[ind] * Texp * 1e4
mask = np.where((sep_s < R_s))
cts1 = np.histogram2d(ev_time[mask], ev_en[mask], bins=[Tbins, EEbins])[0]
mask = sep_b < R_s
cts2 = np.histogram2d(ev_time[mask], ev_en[mask], bins=[Tbins, EEbins])[0]
src = cts1 - cts2
src_err = sqrt(cts1 + cts2)
flux += np.sum(src / AA, axis=1)
flux_err += np.sum(src_err / AA, axis=1)
flux_b += np.sum(cts2 / AA, axis=1)
flux_b_err += np.sum(sqrt(cts2) / AA, axis=1)
mask = EE < Emin
ind_en = len(EE[mask])
Expos += np.histogram(
ev_time,
weights=AA[ind, ind_en] * Texp * 1e4 * np.ones(Nevents) / Nevents,
bins=Tbins,
)[0]
mask = np.where((sep_s < R_s) & (ev_en > Emin) & (ev_en < Emax))
src_cts += np.histogram(ev_time[mask], bins=Tbins)[0]
mask = np.where((sep_b < R_s) & (ev_en > Emin) & (ev_en < Emax))
bkg_cts += np.histogram(ev_time[mask], bins=Tbins)[0]
hdul.close()
print(src_cts)
print(bkg_cts)
print(Expos)
src = src_cts - bkg_cts
src_err = sqrt(src_cts + bkg_cts)
flux = src / (Expos + 1)
flux_err = src_err / (Expos + 1)
flux_b = bkg_cts / (Expos + 1)
flux_b_err = sqrt(bkg_cts) / (Expos + 1)
# flux_err+=np.sum(src_err/AA,axis=1)
# flux_b+=np.sum(cts2/AA,axis=1)
# flux_b_err+=np.sum(sqrt(cts2)/AA,axis=1)
print(flux)

if message == "":
plt.errorbar(
Expand Down Expand Up @@ -193,20 +216,14 @@
)
lc = ODAAstropyTable(Table(data, names=names))

picture = bin_image # http://odahub.io/ontology#ODAPictureProduct
lightcurve_astropy_table = lc # http://odahub.io/ontology#ODAAstropyTable
png = bin_image # http://odahub.io/ontology#ODAPictureProduct
table = lc # http://odahub.io/ontology#ODAAstropyTable

# output gathering
_galaxy_meta_data = {}
_oda_outs = []
_oda_outs.append(("out_Lightcurve_picture", "picture_galaxy.output", picture))
_oda_outs.append(
(
"out_Lightcurve_lightcurve_astropy_table",
"lightcurve_astropy_table_galaxy.output",
lightcurve_astropy_table,
)
)
_oda_outs.append(("out_Lightcurve_png", "png_galaxy.output", png))
_oda_outs.append(("out_Lightcurve_table", "table_galaxy.output", table))

for _outn, _outfn, _outv in _oda_outs:
_galaxy_outfile_name = os.path.join(_galaxy_wd, _outfn)
Expand Down
Loading
Loading