Efficient way to query multiple positions/images; RSP TAP or Butler?

Support
, , ,
#1

Hi folks, as part of the HUN-KON-S1 contribution for performing forced photometry on Solar System moving objects over (potentially years) of past images along the arc, we need to intersect large numbers of positions from an ephemeris with potential Rubin images (using DP0.2 as the proxy currently). We are currently doing this via the RSP TAP to query ivoa.ObsCore and the queries look like this (time/position windows truncated):

SELECT
lsst_visit,
lsst_detector,
lsst_ccdvisitid,
lsst_band,
s_ra,
s_dec,
t_min,
t_max,
s_region
FROM ivoa.ObsCore
WHERE calib_level = 2
AND (lsst_band = ‘r’ OR lsst_band = ‘g’)
AND t_max >= 59802.415775
AND t_min <= 59863.333979
AND (
(t_min >= 59802.415775
AND t_max <= 59802.416122
AND CONTAINS(POINT(‘ICRS’, 58.138352, -41.696853), s_region) = 1)
OR
(t_min >= 59803.360786
AND t_max <= 59803.361133
AND CONTAINS(POINT(‘ICRS’, 58.458592, -41.695724), s_region) = 1)
OR
(t_min >= 59803.362719
AND t_max <= 59803.363066
AND CONTAINS(POINT(‘ICRS’, 58.459217, -41.695733), s_region) = 1)
OR
(t_min >= 59805.314816
AND t_max <= 59805.315163
AND CONTAINS(POINT(‘ICRS’, 59.102676, -41.698524), s_region) = 1)
)

Performance is very slow or just crashes in the query with large numbers of points needed to resolve a long arc.
We’ve investigated other constraints supported by ADQL such as INTERSECTS with POLYGON shapes (since moving objects are often traveling along relatively long but thin paths) based on documentation for the Gaia TAP service. However it seems like support for other than POINTs isn’t available in RSP:

File /opt/lsst/software/stack/conda/envs/lsst-scipipe-10.0.0/lib/python3.12/site-packages/pyvo/dal/tap.py:996, in AsyncTAPJob.raise_if_error(self)
994 msg = self._job.errorsummary.message.content
995 msg = msg or “”
→ 996 raise DALQueryError("Query Error: " + msg, self.url)

DALQueryError: Query Error: UnsupportedOperationException:53 is not a supported coordinate system.

Are there are (documented…?) alternatives using the Butler we should use instead or stick with RSP TAP and find workarounds, better queries or get additional support added (I appreciate “scope creep” is the last thing wanted right now…). Being able to investigate past behavior of newly discovered solar system objects or newly detected outbursts from alerts, is going to be critical for a lot of the time-domain/alert-based solar system science the SSSC wants to do with Rubin.

#2

Hi Tim. We’ll get you a more fulsome reply next week but just to reassure you this is not scope creep, INTERSECT support in Qserv is on our roadmap: Rubin Science Platform Status and Roadmap — Rubin Science Platform

#3

Hi Tim, I noticed you asked about alternatives using the Butler. Since the objective is to do forced photometry, and the butler is optimized for image analysis with the LSST Science Pipelines (like forced photometry), I think a butler-based solution make sense here.

Butler queries do accept statements like detector_region OVERLAPS user_region. It is possible to establish a polygon and a timespan, and query the butler for all visit images (calexps, in DP0.2) in that timespan. Here’s an example that returns 70 DP0.2 visit images (calexps) based on the example constraints in your post.

from lsst.daf.butler import Butler, Timespan
import lsst.sphgeom as sphgeom
from astropy.time import Time

butler = Butler('dp02', collections='2.2i/runs/DP0.2')

polygon = """POLYGON 58.138352 -41.696853
                     58.458592 -41.695724
                     58.459217 -41.695733
                     59.102676 -41.698524"""
region = sphgeom.Region.from_ivoa_pos(polygon)

time1 = Time(59802.415775, format="mjd", scale="tai")
time2 = Time(59863.333979, format="mjd", scale="tai")
timespan = Timespan(time1, time2)
del time1, time2

dataset_refs = butler.query_datasets("calexp",
                                     where="band.name IN ('g', 'r') AND \
                                     visit.timespan OVERLAPS timespan AND \
                                     visit_detector_region.region OVERLAPS region",
                                     bind={"timespan": timespan,
                                           "region": region})
print(len(dataset_refs))

If this example is the kind of thing you were looking for, could you mark this post as the solution? If it wasn’t quite what you needed, we can keep this topic thread open for further discussion.

#4

Hi @talister, can you confirm if @MelissaGraham’s response provided what you were looking for? If so, could you please mark it as the solution? If you still need more information, please provide further clarification and we can keep this thread open for additional discussion.

#5

Sorry @sgreenstreet @MelissaGraham I’m slightly fronting this question for our Hungarian colleagues… I’ve pinged them about @MelissaGraham’s reply and asked if it helped but haven’t heard back yet - in theory we have a stand-up meeting tomorrow, so possibly more news after that…

1 Like
#6

Hi folks, Laszlo (the primary developer on HUN-KON-S1) reports that the new code fragment is working well for him and speeds things up considerably vs TAP, but he has been unable to find any additional details on building and using butler queries beyond the DP02_04 tutorial. Is this something that there could potentially be more documentation for in the future once the First Look and DP1 release excitements calm down ?

#7

See e.g., Querying datasets — LSST Science Pipelines

#8

Thanks Tim. I’m going to mark my reply as the solution but as always, open up a new ticket anytime for new or continuing issues.

Yes, with DP1 we have a few more butler tutorials. It would be interesting to hear more about the types of queries needed and what demonstrations would be most useful, too.

1 Like