Query SOAR FoV Metadata for L2 Remote sensing instruments EUI, PHI and SPICE

[JonCook]

Describe the feature

Recently the solar orbiter archive now has FoV metadata information available for all files via TAP for several remote sensing instruments and sensors:

• EUI L2 HRI174 and HRI1216
• PHI L2 HRT
• SPICE L2

Detailed information about the columns including descriptions is available via the TAP tables page published at:
https://soar.esac.esa.int/soar-sl-tap/tap/tables

• v_spice_fov
• v_eui_hri_fov
• v_phi_hrt_fov

Using standard metadata requests it is possible to retrieve FoV information about a given file or set of files:
For example for a given file:

• http://soar.esac.esa.int/soar-sl-tap/tap/sync?REQUEST=doQuery&LANG=ADQL&FORMAT=json&QUERY=SELECT+*+FROM+v_eui_hri_fov+WHERE+filename=%27solo_L2_eui-hrilya1216-image_20230424T232420166_V01.fits%27

But of course it is also possible to query on any column in the table. There are a few indexes in the columns

With this information it makes it 'quite easy' to make beautiful plots without the knowing the details and internal metadata of each instruments FITS files like:

As well as earth perspective:

And positional:

It could be nice to expose some of this functionality from the soar-sunpy package and also create some nice examples in the galleries

Proposed solution

Consider exposing some of this functionality via the soar-sunpy package and even creating a nice gallery

[hayesla]

This is great @JonCook ! thanks for advertising that these tables are now available.

[nabobalis]

We should make this a gallery example. @NucleonGodX can we add that to the gsoc list?

[JonCook]

If it helps @nabobalis - I have some code which makes the first plot using TAP requests. Perhaps it could be useful to get started. I'm sure it can be cleaned up a lot

[nabobalis]

That would be great! Thank you @JonCook

[NucleonGodX]

We should make this a gallery example. @NucleonGodX can we add that to the gsoc list?

Alright, lets do that.

[JonCook]

That would be great! Thank you @JonCook

Here it is, it is a bit 'rough', but this is the python to produce the first image in the example using TAP requests.. It also assumes you have the FITS files in question downloaded.

#!/usr/bin/env python
# coding: utf-8

# In[2]:


#Set up our dependencies

import sunpy.map
from sunpy.coordinates import frames, get_earth, get_horizons_coord, get_body_heliographic_stonyhurst
from sunpy.net import Fido, attrs as a 
import sunpy_soar
import sunraster.instr.spice
from astropy import units as u 
#from astroquery.utils.tap.core import TapPlus

from astropy.io import fits
from astropy.coordinates import SkyCoord, SkyOffsetFrame
from astropy.time import Time

import pyvo as vo

import matplotlib.pyplot as plt
import math

from dateutil import parser as date_parser
import warnings

warnings.filterwarnings('ignore')

eui_fsi_filename = '/Users/jcook/Downloads/fov-test-files-notebook/solo_L2_eui-fsi174-image_20221024T190050195_V01.fits'
eui_hri_filename = '/Users/jcook/Downloads/fov-test-files-notebook/solo_L2_eui-hrieuv174-image_20221024T190000171_V01.fits'
phi_hrt_filename = '/Users/jcook/Downloads/fov-test-files-notebook/solo_L2_phi-hrt-blos_20221024T191503_V01.fits'
spice_filename   = '/Users/jcook/Downloads/fov-test-files-notebook/solo_L2_spice-n-ras_20221024T190134_V05_150995395-058.fits'

# Lets get some stuff from TAP
service = vo.dal.TAPService("https://soar.esac.esa.int/soar-sl-tap/tap")
fits_item_resulset = service.search("SELECT * FROM soar.v_eui_sc_fits WHERE filename='solo_L2_eui-fsi174-image_20221024T190050195_V01.fits'")
#fits_item_resulset = service.search("SELECT * FROM soar.v_eui_sc_fits")

#print(fits_item_resulset)
print(fits_item_resulset.table.columns)

first_item = fits_item_resulset[0]
fsi_obs_date = first_item["date_average"]
print(first_item["soop_name"])
print(type(first_item["soop_name"]))
print(type(first_item["wavemax"]))

fov_eui_resultset = service.search("SELECT * FROM soar.v_eui_hri_fov WHERE filename='solo_L2_eui-hrieuv174-image_20221024T190000171_V01.fits'")
fov_phi_resultset = service.search("SELECT * FROM soar.v_phi_hrt_fov WHERE filename='solo_L2_phi-hrt-blos_20221024T191503_V01.fits'")
fov_spice_resultset = service.search("SELECT * FROM soar.v_spice_fov WHERE filename='solo_L2_spice-n-ras_20221024T190134_V05_150995395-058.fits'")


# Make a query to get FoV information for a given file
service = vo.dal.TAPService("https://soar.esac.esa.int/soar-sl-tap/tap")
fov_eui_resultset = service.search("SELECT * FROM soar.v_eui_hri_fov WHERE filename='solo_L2_eui-hrieuv174-image_20221024T190000171_V01.fits'")

# Print out the DALResultsTable object (it includes the correct type already)
#print(fov_eui_resultset)
# Print out the TableColumns
#print(fov_eui_resultset.table.columns)

# Print out the soop name and some column types
row = fov_eui_resultset[0]
print(row["soop_name"])
print(type(row["fov_earth_bot_left_arcsec_tx"]))
print(row["fov_earth_bot_left_arcsec_tx"])

eui_fsi_map = sunpy.map.Map(eui_fsi_filename)
eui_hri_map = sunpy.map.Map(eui_hri_filename)
phi_hrt_map = sunpy.map.Map(phi_hrt_filename)

# Generation of SunPy Map for a zoomed version of the FSI image
m_fsi_zoom = eui_fsi_map.submap(SkyCoord(-3000,-3000, unit='arcsec', frame=eui_fsi_map.coordinate_frame),\
                                     top_right=SkyCoord(3000, 3000, unit='arcsec', \
                                     frame=eui_fsi_map.coordinate_frame))

fig = plt.figure(figsize=(20, 10))
ax2 = fig.add_subplot(1, 3, 2, projection=m_fsi_zoom)
## ORBITAL plot
#fsi_obs_date=eui_fsi_info[1].header['DATE-AVG']
m_fsi_zoom.plot(axes=ax2, title='Solar Orbiter perspective (EUI 174 $\mathrm{\AA}$)\n''%s' % fsi_obs_date)

original_bl = SkyCoord(
    fov_eui_resultset["fov_solo_bot_left_arcsec_tx"]*u.arcsec, fov_eui_resultset["fov_solo_bot_left_arcsec_ty"]*u.arcsec, frame=eui_hri_map.coordinate_frame)
original_tr = SkyCoord(
    fov_eui_resultset["fov_solo_top_right_arcsec_tx"]*u.arcsec, fov_eui_resultset["fov_solo_top_right_arcsec_ty"]*u.arcsec, frame=eui_hri_map.coordinate_frame)

rotated_frame = SkyOffsetFrame(origin=eui_hri_map.reference_coordinate,
                                rotation=-fov_eui_resultset["crota"]*u.deg,
                                observer=eui_hri_map.observer_coordinate)

rotated_bl = original_bl.transform_to(rotated_frame)
rotated_tr = original_tr.transform_to(rotated_frame)


# Draw the quadrangle
m_fsi_zoom.draw_quadrangle(
    rotated_bl[0],
    top_right=rotated_tr[0],
    label='EUI HRI',
    edgecolor='C0',
    lw=2)

original_bl = SkyCoord(
    fov_phi_resultset["fov_solo_bot_left_arcsec_tx"]*u.arcsec, fov_phi_resultset["fov_solo_bot_left_arcsec_ty"]*u.arcsec, frame=phi_hrt_map.coordinate_frame)
original_tr = SkyCoord(
    fov_phi_resultset["fov_solo_top_right_arcsec_tx"]*u.arcsec, fov_phi_resultset["fov_solo_top_right_arcsec_ty"]*u.arcsec, frame=phi_hrt_map.coordinate_frame)
rotated_frame = SkyOffsetFrame(origin=phi_hrt_map.reference_coordinate,
                                rotation=-fov_phi_resultset["crota"]*u.deg,
                                observer=phi_hrt_map.observer_coordinate)
rotated_bl = original_bl.transform_to(rotated_frame)
rotated_tr = original_tr.transform_to(rotated_frame)
m_fsi_zoom.draw_quadrangle(
    rotated_bl[0],
    top_right=rotated_tr[0],
    label='PHI HRT',
    edgecolor='C1',
    lw=2)

raster_spice = sunraster.instr.spice.read_spice_l2_fits(spice_filename)
spice_window = raster_spice['O III 703 / Mg IX 706 (Merged)']
# The spice window has data in t, lambda, y, x
spice_map = sunpy.map.Map(spice_window[0, 20, :, :].data, \
                          spice_window[0, 20, :, :].meta)

original_bl = SkyCoord(
    fov_spice_resultset["fov_solo_bot_left_arcsec_tx"]*u.arcsec, fov_spice_resultset["fov_solo_bot_left_arcsec_ty"]*u.arcsec, frame=spice_map.coordinate_frame)
original_tr = SkyCoord(
    fov_spice_resultset["fov_solo_top_right_arcsec_tx"]*u.arcsec, fov_spice_resultset["fov_solo_top_right_arcsec_ty"]*u.arcsec, frame=spice_map.coordinate_frame)
rotated_frame = SkyOffsetFrame(origin=spice_map.reference_coordinate,
                                rotation=-fov_spice_resultset["crota"]*u.deg,
                                observer=spice_map.observer_coordinate)
rotated_bl = original_bl.transform_to(rotated_frame)
rotated_tr = original_tr.transform_to(rotated_frame)
m_fsi_zoom.draw_quadrangle(
    rotated_bl[0],
    top_right=rotated_tr[0],
    label='SPICE',
    edgecolor='C2',
    lw=2)

ax2.legend()

Thanks here to @hayesla and others. I just changed a bit here where we get the values from. We also have the code for generating the earth view and positional images, but I did not update it yet to work with TAP

Thanks

[JonCook]

Hi @nabobalis and @NucleonGodX - just wondered if work continues here or it got frozen? Thanks

[nabobalis]

Hi @JonCook, the GSoC project has ended so I am not sure any more work will continue on this in the short term.

We did this issue: #143 (comment) but we were not able to come up with a solution.