fov table joins added

examples/field_of_view_plot.py

@@ -0,0 +1,84 @@
+"""
+========================================
+Retrieve and Plot Field of View on a map
+========================================
+
+This example demonstrates how to fetch Solar Orbiter's Field of View (FOV) for each instrument values
+and Solar Orbiter data using `sunpy.net.Fido` and plot them on a `sunpy.map.Map`.
+"""
+
+import astropy.units as u
+import matplotlib.pyplot as plt
+import numpy as np
+import sunpy.map
+import sunpy.net.attrs as a
+from astropy.coordinates import SkyCoord
+from sunpy.coordinates import frames
+from sunpy.net import Fido
+
+#####################################################
+# Importing sunpy_soar registers the client with sunpy
+import sunpy_soar  # NOQA: F401
+
+#####################################################
+# We'll begin by constructing a search query for the Field of View (FOV).
+# In this example, we'll use the instrument EUI and set the perspective to that of Earth's field of view.
+
+instrument = a.Instrument("EUI")
+time = a.Time("2022-10-13 12:06:00", "2022-10-13 12:06:10")
+level = a.Level(2)
+detector = a.Detector("HRI_EUV")
+product = a.soar.Product("eui-hrieuv174-image")
+fov = a.soar.FOV("earth")
+
+#####################################################
+# Now we will do the search.
+
+result = Fido.search(instrument & time & level & detector & product & fov)
+
+#####################################################
+# To plot the FOV, we need a map to overlay them on to.
+# For this we will create a blank map as the base.
+
+data = np.full((10, 10), np.nan)
+skycoord = SkyCoord(
+    0 * u.arcsec, 0 * u.arcsec, obstime="2022-10-13 12:06:00", observer="earth", frame=frames.Helioprojective
+)
+header = sunpy.map.make_fitswcs_header(data, skycoord, scale=[220, 220] * u.arcsec / u.pixel)
+blank_map = sunpy.map.Map(data, header)
+
+#####################################################
+# Now we need to  extract the FOV coordinates from search results.
+
+fov_earth_bot_left_tx = result[0]["fov_earth_left_arcsec_tx"][0] * u.arcsec
+fov_earth_bot_left_ty = result[0]["fov_earth_left_arcsec_ty"][0] * u.arcsec
+fov_earth_top_right_tx = result[0]["fov_earth_right_arcsec_tx"][0] * u.arcsec
+fov_earth_top_right_ty = result[0]["fov_earth_right_arcsec_ty"][0] * u.arcsec
+
+# Ensure correct ordering of coordinates
+bottom_left_tx = min(fov_earth_bot_left_tx, fov_earth_top_right_tx)
+bottom_left_ty = min(fov_earth_bot_left_ty, fov_earth_top_right_ty)
+top_right_tx = max(fov_earth_bot_left_tx, fov_earth_top_right_tx)
+top_right_ty = max(fov_earth_bot_left_ty, fov_earth_top_right_ty)
+
+#####################################################
+# To plot the corners of the corners of the FOVs, we need turn them into a `~astropy.coordinates.SkyCoord`.
+
+earth_fov_bottom_left = SkyCoord(bottom_left_tx, bottom_left_ty, frame=blank_map.coordinate_frame)
+earth_fov_top_right = SkyCoord(top_right_tx, top_right_ty, frame=blank_map.coordinate_frame)
+
+# Plot the blank map and the FOV.
+fig = plt.figure()
+ax = fig.add_subplot(projection=blank_map)
+blank_map.plot(axes=ax)
+blank_map.draw_limb(axes=ax, color="k")
+blank_map.draw_grid(axes=ax, color="k")
+
+# Draw the FOV as a quadrangle
+blank_map.draw_quadrangle(earth_fov_bottom_left, top_right=earth_fov_top_right, axes=ax, edgecolor="blue")
+
+# Set title and show legend
+ax.set_title("Fields of View on Blank Map")
+ax.legend()
+
+plt.show()

pyproject.toml

@@ -19,6 +19,7 @@ authors = [
 dependencies = [
   "astropy>=5.3.0",
   "sunpy[net]>=6.0.0",
+  "sunpy[image]",
   "requests>=2.28.0",
 ]
 dynamic = ["version"]

sunpy_soar/attrs.py

@@ -19,6 +19,14 @@ def __init__(self, value):
         self.value = value.lower()
 
 
+class FOV(SimpleAttr):
+    """
+    The Field of View type to fetch coordinates for.
+
+    Allowed values are "solar" and "earth".
+    """
+
+
 class SOOP(SimpleAttr):
     """
     The SOOP name to search for.
@@ -141,3 +149,15 @@ def _(wlk, attr, params):  # NOQA: ARG001
     wavemin = attr.min.value
     wavemax = attr.max.value
     params.append(f"Wavemin='{wavemin}'+AND+Wavemax='{wavemax}'")
+
+
+@walker.add_applier(FOV)
+def _(wlk, attr, params):  # NOQA: ARG001
+    allowed_fovs = ["solar", "earth"]
+    if attr.value not in allowed_fovs:
+        warnings.warn(
+            f"FOV not in list of allowed FOVs for SOAR: {allowed_fovs}",
+            SunpyUserWarning,
+            stacklevel=2,
+        )
+    params.append(f"FOV ='{attr.value}'")

sunpy_soar/client.py

@@ -12,7 +12,7 @@
 from sunpy.net.base_client import BaseClient, QueryResponseTable
 from sunpy.time import parse_time
 
-from sunpy_soar.attrs import SOOP, Product, walker
+from sunpy_soar.attrs import FOV, SOOP, Product, walker
 
 __all__ = ["SOARClient"]
 
@@ -34,10 +34,88 @@ def search(self, *query, **kwargs):  # NOQA: ARG002
         table = astropy.table.vstack(results)
         qrt = QueryResponseTable(table, client=self)
         qrt["Filesize"] = (qrt["Filesize"] * u.byte).to(u.Mbyte).round(3)
-        qrt.hide_keys = ["Data item ID", "Filename"]
+        qrt.hide_keys = [
+            "Data item ID",
+            "Filename",
+        ]
         return qrt
 
-    def add_join_to_query(query: list[str], data_table: str, instrument_table: str):
+    @staticmethod
+    def determine_fov_table(query):
+        """
+        Determine the value of 'm' based on the FOV parameter in the query.
+
+        Parameters
+        ----------
+        query : list[str]
+            List of query items.
+
+        Returns
+        -------
+        str or None
+            The value of 'm' ("earth", "solo", or None).
+        """
+        for q in query:
+            match = re.match(r"FOV\s*=\s*'([^']+)'", q)
+            if match:
+                fov_value = match.group(1).lower()
+                if fov_value == "earth":
+                    return "earth"
+                if fov_value == "solar":
+                    return "solo"
+        return None
+
+    def fov_join(self, query, instrument_table: str, where_part: str, from_part: str, select_part: str):
+        """
+        Add FoV (Field of View) join to the query if applicable based on the
+        instrument.
+
+        Parameters
+        ----------
+        instrument_table : str
+            Name of the instrument table.
+        where_part : str
+            The WHERE part of the ADQL query.
+        from_part : str
+            The FROM part of the ADQL query.
+        select_part : str
+            The SELECT part of the ADQL query.
+
+        Returns
+        -------
+        tuple[str, str, str]
+            Updated WHERE, FROM, and SELECT parts of the query.
+        """
+        m = self.determine_fov_table(query)
+        if not m:
+            return where_part, from_part, select_part
+
+        join_tables = {
+            "eui": (
+                "v_eui_hri_fov",
+                f"fov_{m}_bot_left_arcsec_ty, h3.fov_{m}_bot_left_arcsec_tx, "
+                f"h3.fov_{m}_top_right_arcsec_ty, h3.fov_{m}_top_right_arcsec_tx",
+            ),
+            "spi": (
+                "v_spice_fov",
+                f"fov_{m}_bot_left_arcsec_ty, h3.fov_{m}_bot_left_arcsec_tx, "
+                f"h3.fov_{m}_top_right_arcsec_ty, h3.fov_{m}_top_right_arcsec_tx",
+            ),
+            "phi": (
+                "v_phi_hrt_fov",
+                f"fov_{m}_bot_left_arcsec_ty, h3.fov_{m}_bot_left_arcsec_tx, "
+                f"h3.fov_{m}_top_right_arcsec_ty, h3.fov_{m}_top_right_arcsec_tx",
+            ),
+        }
+        for instrument, (table, fields) in join_tables.items():
+            if instrument in instrument_table:
+                from_part += f" LEFT JOIN {table} AS h3 ON h2.filename = h3.filename"
+                select_part += f", {fields}"
+                break
+
+        return where_part, from_part, select_part
+
+    def add_join_to_query(self, query: list[str], data_table: str, instrument_table: str):
         """
         Construct the WHERE, FROM, and SELECT parts of the ADQL query.
 
@@ -60,37 +138,39 @@ def add_join_to_query(query: list[str], data_table: str, instrument_table: str):
         wavemin_pattern = re.compile(r"Wavemin='(\d+\.\d+)'")
         wavemax_pattern = re.compile(r"Wavemax='(\d+\.\d+)'")
         for parameter in query:
+            if parameter.startswith("FOV"):
+                continue
             wavemin_match = wavemin_pattern.search(parameter)
             wavemax_match = wavemax_pattern.search(parameter)
-            # If the wavemin and wavemax are same that means only one wavelength is given in query.
+            # If the wavemin and wavemax are the same, that means only one wavelength is given in the query.
             if wavemin_match and wavemax_match and float(wavemin_match.group(1)) == float(wavemax_match.group(1)):
-                # For PHI and SPICE, we can specify wavemin and wavemax in the query and get the results.
-                # For PHI we have wavelength data in both angstrom and nanometer without it being mentioned in the SOAR.
-                # For SPICE we get data in form of wavemin/wavemax columns, but only for the first spectral window.
-                # To make sure this data is not misleading to the user we do not return any values for PHI AND SPICE.
                 parameter = f"Wavelength='{wavemin_match.group(1)}'"
             elif wavemin_match and wavemax_match:
                 parameter = f"Wavemin='{wavemin_match.group(1)}'+AND+h2.Wavemax='{wavemax_match.group(1)}'"
-            prefix = "h1." if not parameter.startswith("Detector") and not parameter.startswith("Wave") else "h2."
+            prefix = "h2." if parameter.startswith(("Detector", "Wave", "Observation")) else "h1."
             if parameter.startswith("begin_time"):
                 time_list = parameter.split("+AND+")
                 final_query += f"h1.{time_list[0]}+AND+h1.{time_list[1]}+AND+"
-                # As there are no dimensions in STIX, the dimension index need not be included in the query for STIX.
                 if "stx" not in instrument_table:
-                    # To avoid duplicate rows in the output table, the dimension index is set to 1.
                     final_query += "h2.dimension_index='1'+AND+"
             else:
                 final_query += f"{prefix}{parameter}+AND+"
 
         where_part = final_query[:-5]
         from_part = f"{data_table} AS h1"
         select_part = (
-            "h1.instrument, h1.descriptor, h1.level, h1.begin_time, h1.end_time, "
-            "h1.data_item_id, h1.filesize, h1.filename, h1.soop_name"
+            "h1.instrument, h1.descriptor, h1.level, h1.begin_time, h1.end_time, h1.data_item_id, "
+            "h1.filesize, h1.filename, h1.soop_name, h2.wavelength, h2.detector, h2.dimension_index"
         )
-        if instrument_table:
-            from_part += f" JOIN {instrument_table} AS h2 USING (data_item_oid)"
-            select_part += ", h2.detector, h2.wavelength, h2.dimension_index"
+        # Add the second join always
+        from_part += f" JOIN {instrument_table} AS h2 ON h1.data_item_oid = h2.data_item_oid"
+
+        # Add the third join conditionally based on the instrument and other conditions
+        if any(q.startswith("FOV") for q in query):
+            where_part, from_part, select_part = self.fov_join(
+                query, instrument_table, where_part, from_part, select_part
+            )
+
         return where_part, from_part, select_part
 
     @staticmethod
@@ -137,7 +217,7 @@ def _construct_payload(query):
 
         # Need to establish join for remote sensing instruments as they have instrument tables in SOAR.
         if instrument_name in ["EUI", "MET", "SPI", "PHI", "SHI"]:
-            where_part, from_part, select_part = SOARClient.add_join_to_query(query, data_table, instrument_table)
+            where_part, from_part, select_part = SOARClient().add_join_to_query(query, data_table, instrument_table)
         else:
             from_part = data_table
             select_part = "*"
@@ -167,7 +247,6 @@ def _do_search(query):
         payload = SOARClient._construct_payload(query)
         # Need to force requests to not form-encode the parameters
         payload = "&".join([f"{key}={val}" for key, val in payload.items()])
-        # Get request info
         r = requests.get(f"{tap_endpoint}/sync", params=payload)
         log.debug(f"Sent query: {r.url}")
         r.raise_for_status()
@@ -189,23 +268,39 @@ def _do_search(query):
             info["begin_time"] = parse_time(info["begin_time"]).iso
             info["end_time"] = parse_time(info["end_time"]).iso
 
-        result_table = astropy.table.QTable(
-            {
-                "Instrument": info["instrument"],
-                "Data product": info["descriptor"],
-                "Level": info["level"],
-                "Start time": info["begin_time"],
-                "End time": info["end_time"],
-                "Data item ID": info["data_item_id"],
-                "Filename": info["filename"],
-                "Filesize": info["filesize"],
-                "SOOP Name": info["soop_name"],
-            },
-        )
-        if "detector" in info:
-            result_table["Detector"] = info["detector"]
-        if "wavelength" in info:
-            result_table["Wavelength"] = info["wavelength"]
+        m = SOARClient.determine_fov_table(query)
+        contains_fov = any(q.lower().startswith("fov") for q in query)
+        if not contains_fov:
+            result_table = astropy.table.QTable(
+                {
+                    "Instrument": info["instrument"],
+                    "Data product": info["descriptor"],
+                    "Level": info["level"],
+                    "Start time": info["begin_time"],
+                    "End time": info["end_time"],
+                    "Data item ID": info["data_item_id"],
+                    "Filename": info["filename"],
+                    "Filesize": info["filesize"],
+                    "SOOP Name": info["soop_name"],
+                },
+            )
+            if "wavelength" in info:
+                result_table["Wavelength"] = info["wavelength"]
+            if "detector" in info:
+                result_table["Detector"] = info["detector"]
+        else:
+            result_table = astropy.table.QTable(
+                {
+                    "Instrument": info["instrument"],
+                    "Start time": info["begin_time"],
+                    "End time": info["end_time"],
+                    "Filesize": info["filesize"],
+                    f"fov_{m}_left_arcsec_ty": info[f"fov_{m}_bot_left_arcsec_ty"],
+                    f"fov_{m}_left_arcsec_tx": info[f"fov_{m}_bot_left_arcsec_tx"],
+                    f"fov_{m}_right_arcsec_ty": info[f"fov_{m}_top_right_arcsec_ty"],
+                    f"fov_{m}_right_arcsec_tx": info[f"fov_{m}_top_right_arcsec_tx"],
+                },
+            )
         result_table.sort("Start time")
         return result_table
 
@@ -253,7 +348,7 @@ def _can_handle_query(cls, *query):
             True if this client can handle the given query.
         """
         required = {a.Time}
-        optional = {a.Instrument, a.Detector, a.Wavelength, a.Level, a.Provider, Product, SOOP}
+        optional = {a.Instrument, a.Detector, a.Wavelength, a.Level, a.Provider, Product, SOOP, FOV}
         if not cls.check_attr_types_in_query(query, required, optional):
             return False
         # check to make sure the instrument attr passed is one provided by the SOAR.

sunpy_soar/tests/test_sunpy_soar.py

@@ -234,6 +234,17 @@ def test_wavelength_range():
         assert all(table["Wavelength"] == 174)
 
 
+def test_fov_search():
+    instrument = a.Instrument("SPICE")
+    time = a.Time("2022-10-12 16:00", "2022-10-12 16:20")
+    level = a.Level(2)
+    product = a.soar.Product("spice-n-sit")
+    Fov = a.soar.FOV("solar")
+
+    res = Fido.search(instrument & time & level & product & Fov)
+    assert res[0]["fov_solo_left_arcsec_tx"][0] == 16.800173739034108
+
+
 def test_join_science_query():
     result = SOARClient._construct_payload(  # NOQA: SLF001
         [
@@ -246,8 +257,8 @@ def test_join_science_query():
 
     assert result["QUERY"] == (
         "SELECT+h1.instrument, h1.descriptor, h1.level, h1.begin_time, h1.end_time, "
-        "h1.data_item_id, h1.filesize, h1.filename, h1.soop_name, h2.detector, h2.wavelength, "
-        "h2.dimension_index+FROM+v_sc_data_item AS h1 JOIN v_eui_sc_fits AS h2 USING (data_item_oid)"
+        "h1.data_item_id, h1.filesize, h1.filename, h1.soop_name, h2.wavelength, h2.detector, "
+        "h2.dimension_index+FROM+v_sc_data_item AS h1 JOIN v_eui_sc_fits AS h2 ON h1.data_item_oid = h2.data_item_oid"
         "+WHERE+h1.instrument='EUI'+AND+h1.begin_time>='2021-02-01+00:00:00'+AND+h1.begin_time<='2021-02-02+00:00:00'"
         "+AND+h2.dimension_index='1'+AND+h1.level='L1'+AND+h1.descriptor='eui-fsi174-image'"
     )
@@ -265,8 +276,8 @@ def test_join_low_latency_query():
 
     assert result["QUERY"] == (
         "SELECT+h1.instrument, h1.descriptor, h1.level, h1.begin_time, h1.end_time, "
-        "h1.data_item_id, h1.filesize, h1.filename, h1.soop_name, h2.detector, h2.wavelength, "
-        "h2.dimension_index+FROM+v_ll_data_item AS h1 JOIN v_eui_ll_fits AS h2 USING (data_item_oid)"
+        "h1.data_item_id, h1.filesize, h1.filename, h1.soop_name, h2.wavelength, h2.detector, "
+        "h2.dimension_index+FROM+v_ll_data_item AS h1 JOIN v_eui_ll_fits AS h2 ON h1.data_item_oid = h2.data_item_oid"
         "+WHERE+h1.instrument='EUI'+AND+h1.begin_time>='2021-02-01+00:00:00'+AND+h1.begin_time<='2021-02-02+00:00:00'"
         "+AND+h2.dimension_index='1'+AND+h1.level='LL01'+AND+h1.descriptor='eui-fsi174-image'"
     )